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Background 


It has long been recognized that subtle changes in the orbital and rotational parameters of the Earth 
directly modify the spatial and temporal patterns of incident solar radiation, and thus indirectly, but still 
quite significantly, influence climate. Though the time scales of climatic forcing via orbital and rotational 
change are far removed from the normal lime scales of human awareness and public policy decisions, they 
arc nonetheless important The importance of astronomically forced climatic change stems from two 
considerations: the overall magnitude of the effect, and the potential for calibration of climate system 
response on other time scales. The dominant climatic signal on time scales from 10 4 years to 10 7 years is 
phase coherent with, and presumably causally linked to orbital and rotational variations. Because the 
astronomical forcing can be computationally reconstructed with considerable accuracy, and the climatic 
response is preserved in a myriad of proxy records, the system response function in this frequency band 
can be quite well constrained. 


Objectives 


The basic scientific goal of the Workshop was to arrive at a better understanding of the interactions 
between the orbital, rotational and climatic variations of the Earth. In order to accomplish that objective, 
the Workshop was structured so as to provide an opportunity for experts in the various sub-disciplines to 
vigorously interact and exchange ideas, in a rather informal setting. 

A number of previous workshops and symposia have been held over the last decade on various aspects of 
orbital, rotational and climatic interactions. However, significant recent developments have taken place in 
a number of the scientific disciplines relevant to understanding the mechanisms associated with 
astronomical forcing and climatic response. As the field matures, discussion can shift from the issues 
associated with "whether" to those associated with "how". Some of these recent developments include: 
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improved models of planetary orbital evolution over 10 to 10 year timespans, 

enhanced understanding of the Earth’s defonnational and rotational response to glacial loading 

better appreciation for the role of core-mantle coupling in precessional dynamics 

improved radiometric dating of key events of the last glacial cycle 

acquisition of additional marine and continental climate records 

implementation of improved statistical techniques for comparing climatic and astronomical time series 
maturation of algorithms for modeling climatic states and processes 

Until quite recently, most of the people actively working on astronomically forced climatic change were 
oceanographers or glaciologists. The simple explanation for this observation is that the best records of 
quasiperiodic climatic change came from marine sediments and glacial ice cores. However, as the field 
develops, one can expect participation from a wider range of disciplines. For example, it seems likely that 
significant progress can be made in understanding the precessional dynamics of the solid Earth by 
examining the paleocliniatic record. Though recent progress in space geodetic techniques has made it 
possible to precisely characterize the rotational response of the Earth on time scales from a few hours to 
the lunar nodal precession period (18.6 years), the details of the Earth's response on the spin precession 
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time scale (26 thousand years) are still poorly known, since the direct observations only span a short part 
of one period. In that case, high resolution spectral studies of long climatic records have already indicated 
subtle effects associated with glacial mass redistribution. 

Because of the very interdisciplinary nature of the subject matter of the Workshop, it is difficult for any 
one person, or any small group of people to be really conversant with the latest results, techniques and 
opinions across the entire spectrum of relevant topics. However, by assembling a diverse and reasonably 
large group of people with complementary interests, it was anticipated that each of the participants would 
gain a better appreciation for those areas outside their immediate expertise. It was hoped that each 
participant would come away from the meeting with an improved understanding of symplectic integration 
schemes, benthic foraminifera! taphonomy, isotopic fractionation dynamics, stochastic resonances, visco- 
elastic normal modes, Hopf bifurcations, quasi-geostrophic potential vorticity, and many of the other 
topics represented. 

However, the objectives of the workshop extended beyond educating the individual participants. It was 
also anticipated that collectively, we would take a critical look at what is now known in various specialty 
areas, and also take a broad look at how well all the little pieces fit together. This kind of examination is 
helpful in identifying promising new research directions which will eventually fill in the gaps in our 
present understanding. 


Organization 


On 9-11 July 1991, an international workshop on Orbital, Rotational and Climatic Interactions was held 
in Olin Hall, which houses the Department of Earth and Planetary Sciences, at the Johns Hopkins 
University, in Baltimore, Maryland. The workshop was hosted by the NASA Goddard Space Flight 
Center, and was convened by Bruce G. Bills, Head of the Geodynamics Branch at GSFC. The agenda 
consisted of five half-day sessions, in which each participant was given an opportunity to present a 
summary of their recent work, with ample lime for discussion. The afternoon of the third day was devoted 
to informal discussion of directions for future research. 


Acknowledgments 

All workshop participants are thanked for the contributions to this report. Many others also provided 
valuable input. Financial support was provided by the NASA Biogeochemistiy and Geophysics Program. 



Workshop Participants 


Bruce G. Bills 
Geodynamics Branch 
Goddard Space Flight Center 
Greenbelt, MD 20771 


Gerard Bond 

Lamont-Doherty Geological Observatory' 
Columbia University 
Palisades, NY 10964 


Madeleine Briskin 
Department of Geology 
University' of Cincinnati 
Cincinnati, OH 45221 


Benjamin F. Chao 
Geodynamics Branch 
Goddard Space Flight Center 
Greenbelt, MD 2077 1 


Steven C. Clemens 
Department of Geological Sciences 
Brown University 
Providence, RI 02912 


Thomas J. Crowley 
Applied Research Corporation 
305 Arguello Drive 
College Station, TX 77840 


Alfred G. Fischer 
Department of Geological Sciences 
University of Southern California 
Los Angeles, CA 90089 


L. Danny Harvey 
Department of Geography 
University' of Toronto 
Toronto, Ontario 
M5S 1A1, Canada 

Timothy D. Herbert 
Scripps Institution of Oceanography 
University of California 
La Jolla, CA 92093 


Linda Hinnov 

Department of Earth and Planetary Sciences 
The Johns Hopkins University 
Baltimore, MD 21218 


William T. Hyde 
Department of Oceanography 
Dalhousie University 
Halifax, Nova Scotia 
B3K 5T8, Canada 

Jacques Laskar 

Service de Mechanique Celeste 
Bureau des Longitudes 
77 Avenue Denfert-Rochereau 
F75014 Paris, France 

Mane-France Loutre 

Institut dAstonomie et Geophysique 

Universite Catholique de Louvain 

2, Chemin du Cyclotron 

B1348 Louvain-la-Neuve, Belgium 

Andrew McIntyre 

Lamont-Doherty Geological Observatory 
Columbia University 
Palisades, NY 10964 


Paul E. Olsen 

Lamont-Doherty Geological Observatory 
Columbia University 
Palisades, NY 10964 


Jeffrey Park 

Department of Geology and Geophysics 
Yale University 
New Haven, CT 065 1 1 


W. Richard Peltier 
Department of Physics 
University of Toronto 
Toronto, Ontario 
M5S 1A7, Canada 

David L. Richardson 
Department of Aerospace Engineering 
University of Cincinnati 
Cincinnati, OH 45221 


vii 



Barry Saltzman 

Department of Geology and Geophysics 
Yale University 
New Haven, CT 06511 


David J. Thomson 

Mathematical Sciences Research Center 
AT&T Bell Laboratories 
Murray Hill, NJ 07974 


Kenneth L. Verosub 
Geology Department 
University of California 
Davis, CA 95616 


Patrick Wu 

Department of Geology and Geophysics 
University of Calgary' 

Calgary, Alberta 
T2N 1N4, Canada 


N93-18933 


Geodynamic Contributions to Global Climatic Change 


Bruce G. Bills 


Geodynamics Branch 
NASA Goddard Space Flight Center 
Greenbelt, MD 20771 


Abstract 

Orbital and rotational variations perturb the latitudinal and seasonal pattern of incident 
solar radiation, producing major climatic change on time scales of 10 4 -10 6 years. The orbital 
variations are oblivious to internal structure and processes, but the rotational variations are not. 
The intent of this article is to describe a program of investigation whose objective would be to 
explore and quantify three aspects of orbital, rotational and climatic interactions. An important 
premise of this investigation is the synergism between geodynamics and paleoclimate. Better 
geophysical models of precessional dynamics are needed in order to accurately reconstruct the 
radiative input to climate models. Some of the paleoclimate proxy records contain information 
relevant to solid Earth processes, on time scales which are difficult to constrain otherwise. 

Specific mechanisms which will be addressed include: 

• climatic consequences of deglacial polar motion, and 

• precessional and climatic consequences of 

- glacially induced perturbations in the gravitational oblateness, and 

— partial decoupling of the mantle and core. 

The approach entails constructing theoretical models of the rotational, deformational, ra- 
diative and climatic response of the Earth to known orbital perturbations, and comparing these 
with extensive records of paleoclimate proxy data. Several of the mechanisms of interest may 
participate in previously unrecognized feed-back loops in the climate dynamics system. A new 
algorithm for estimating climatically diagnostic locations and seasons from the paleoclimate 
time series is proposed. 
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INTRODUCTION 


It is widely recognized that mass redistributions associated with climatic change (glaciations) 
are an important source of crustal deformation and geodynamic change. It is much less widely 
appreciated that rates, phases and amplitudes of deformation of the deep interior of the Earth 
can influence climate. The objective of this investigation is to better characterize four aspects of 
this geodynamic contribution to global climatic change. The common theme among them has two 
threads: internal mass redistributions influence the rotational dynamics of the Earth, and changes 
in orbital arid rotational parameters influence the latitudinal and seasonal pattern of insolation. 
Previous attempts to account for astronomically forced climatic change have usually only considered 
extremely simplisitic models for the response of the Earth to external torques and surface loads. 

The focus of the proposed investigation will consist of three specific aspects of interaction between 
internal mass flow and rotational dynamics, and a new algorithm for comparing paleoclimatic proxy 
data with astronomically forced variations in insolation patterns. 

The latitudinal and seasonal pattern of incident solar radiation depends on the eccentricity of the 
Earth’s orbit and the orientation of the spin axis relative to both the orbit normal and the aspidal 
line. Unit vectors s and n characterize the directions of the spin axis and orbit normal, respectively. 
Two angles completely characterize the relative orientation of the spin axis. The obliquity e is simply 
the angle between the orbit normal and the spin axis 

e = cos~ l (n • s) (1) 

The ascending node of the orbit plane on the instantaneous equator plane has an orientation given 
by (s x n), and the longitude of perihelion w is just the angle in the orbit plane from that node to 
perihelion. It is widely appreciated that secular variations in these three parameters (e,e,r^) produce 
major climatic change (Hayset al., 1976; Berger et ah, 1984). In fact, spectral analyses of long, high 
resolution marine sediment isotopic records show significant variance at periods near 100 kyr, 41 kyr 
and 19-23 kyr, which are generally attributed to spectral lines in the radiative forcing fluctuations 
associated with e, c and e sin(ti7), respectively. 

The causes and effects of the orbital changes are quite well understood. Gravitational interactions 
with the other planets cause the shape and orientation of the orbit to change on time scales of 
10 4 -10 6 years. The inclination I and nodal longitude fi determine the orientation of the orbit plane. 
The eccentricity e and perihelic longitude u; determine the shape of the orbit and its orientation 
within the plane. Note that w is measured from an inertially fixed direction, rather than the moving 
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node as is the case for w. The secular evolution of the orbital element pairs (I, ft) and (e,w) can be 
conveniently represented in terms of Poisson series 


p = 

sin(I) sin(Q) 

= E Nj sin(sj t + gj) 

(2) 

Q = 

sin(I) cos(^) 

= E Nj cos(sj t + gj) 


h 

= e sin(u) = 

E Mj sin(rj t + fj) 

(3) 

k 

= e cos(u>) = 

E Mj cos(r } t + fj) 


In the lowest order solution, there are as many frequencies r, and Sj as there are planets. However, 
the frequencies r ; and Sj are characteristic modal frequencies (eigenvalues) of the coupled system 
of oscillators and are not each uniquely associated with a particular planet (Milam, 1988). The 
frequencies r ; - are all positive, indicating that the perihelia advance. In the lowest order solution, 
the apsidal rates are all in the interval (0.667 < r, < 28.221 arcsec/year). The corresponding 
periods are 45.92 kyr to 1.943 Myr. One of the frequencies s is zero, and all the others are negative, 
indicating that the nodes regress. In the lowest order solution, the non-zero nodal rates are all in 
the interval (0.692 < Sj < 26.330 arcsec/year). The corresponding periods are 49.22 kyr to 1.873 
Myr. In higher order solutions, variations in (e,w) become coupled to variations in (I, ft), but the 
solutions can still be cast in terms of Poisson series like Equations (2) and (3). 

Laskar (1988) has recently published a secular variation theory which is complete to fifth order in 
eccentricity and inclination. Agreement between this secular variation model and strictly numerical 
computations (Richardson and Walker, 1989; Quinn et al„ 1991) is much better than for any previous 
analytical model. The inclination and eccentricity series for Earth each contain 80 distinct terms. 
Figures (1) and (2) illustrate the spectra and corresponding histories of variation in orbital inclination 
and eccentricity for the past 2 10 6 years. 

In computing these secular orbital variations, the Earth, Moon and planets can all be treated as 
point masses. No internal structure or processes are relevant to orbital evolution. The physics of 
the process is simple, and well understood, though development of proper mathematical tools to 
represent the long term evolution remains an area of active research (Laskar, 1988, 1990; Quinn et 
ah, 1991). On the other hand, the rotational evolution does depend rather sensitively on various 

aspects of the structure and dynamics of the interior. 

Lunar and solar gravitational torques acting on the oblate figure of the Earth cause the spin axis 
s to precess about the instantaneous orbit normal n. If the Earth is considered to be a rigid body, 
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the evolution of the spin axis orientation is given by 


ds/dt = o(n • s)(s x n) ( 4 ) 

where 

_ 3 (C — A) ^ Gm, , 9 

“ " “ 2Cn E - 3sm ( 7 ')> (5) 

is a scalar rate factor which depends on intrinsic properties of the Earth, such as polar and 
equatorial moments of inertia (C,A) and rotation rate n, and on extrinsic influences, such as masses 
m, orbital inclinations I, and semiminor axes b, of the Moon and Sun. The solar and lunar torques 
together produce a precession of the spin axis of the Earth at a rate of a(n ■ s) = 50.38 arcsec/year 
(Kinoshita, 1977, Williams et al., 1991). 

Once the present spin axis direction s is known and orbital element histories are given via Equations 
(2) and (3), an obliquity history can be constructed from equation ( 4 ) in two different ways. The 
linear perturbation approach (Miskovic, 1931, Sharaf and Boudnikova, 1967; Vernekar, 1972, 1977; 
Ward, 1974; Berger, 1976) involves deriving coefficients of a trigonometric series, similar to Equations 
(2) and (3), which yield the obliquity and longitude of perihelion directly as functions of time. An 
alternative is to apply standard numerical algorithms for solving initial value problems to generate 
a vector time series s(t) and then compute the obliquity and longitude of perihelion directly (Ward, 
1979; Laskar, 1986; Bills, 1990b). Figure (3a) shows the spectrum of obliquity variations, which 
in the linear perturbation model is simply obtained from the inclination spectrum by shifting each 
frequency Sj by the luni-solar precession rate (s = 50.38 arcsec/year) and multiplying each amplitude 
Nj by the spectral admittance 

Fj = s } /(s } + s) (6) 

Figure (3b) shows the numerically integrated obliquity history of Laskar (1986,1988) for the last 
2 Ma, using the rotational theory of Kinoshita (1977). The difference between the numerical and 
linear perturbation solution never exceeds 0.06 degree over that interval. It is clear that the linear 
perturbation solution gives a very adequate representation of the spin precession. 

The influence of orbital and rotational variations on climate is operative through perturbations 
in the latitudinal and seasonal pattern of insolation. The diurnal average intensity of radiation at a 
point is inversely proportional to the squared solar distance and directly proportional to the diurnal 
average rectified solar direction cosine 


F = (°/ r ) 2 (ll u ' U,||) 


(7) 
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where a and r are mean and instantaneous solar distance and u and u, are unit vectors from 
the center of the Earth to the surface point of interest and the sub-solar point, respectively. The 
insolation pattern, as a function of latitude 6 and mean anomaly M, can be readily computed 
once values are specified for the orbital and rotational parameters c, e, and w (Hargreaves, 1895; 
Milankovitch, 1920; Vernekar, 1972, 1977; Ward, 1974). Figure (4) illustrates that pattern for 
the present orbital and rotational configuration. This pattern can also be written in terms of a 
Fourier-Legendre series (Hargreaves, 1895; North and Coakley, 1979; Taylor, 1984; Bills, 1992) 

F(fi, ^ P n {v) Yh ex P( l P M ) F n , p ( c ' e - ro ) ( 8 ) 

where // = cos($), and P n is a Legendre polynomial. The number of terms in the Fourier summa- 
tion required to obtain a good representation of the seasonal pattern is greater in the polar regions 
than in the tropics and mid-latitudes. The primary difficulty in the polar regions is reproducing 
the abrupt change in slope of the insolation curve at times of transition to continual darkness or 
continual light. It is also true that the polar regions place the greatest demands on the Legendre 
summation, since the spatial pattern also has a discontinuous first derivative at the latitude where 
the transition occurs to continual darkness or light. 

A significant fraction of the recent work on comparing paleoclimate proxy records to astronomically 
forced insolation changes has been based on the published insolation curves of Berger (1978, 1991). 
Common practice is to compare computed variations in some particular aspect of the seasonal and 
latitudinal insolation pattern (July insolation at 65° N is a particularly frequent choice) with an 
observational record of some climatic indicator (<5 18 0 variations versus age (depth) in a marine 
sediment core, for example). Comparisons of this sort enable estimates of amplitude, phase and 
coherence of climatic response to radiative forcing, and have unequivocally demonstrated that orbital 
and rotational variations are a dominant cause of climatic change on time scales of 10 -10 years. 

Despite obvious successes, several problems remain in the general methodology. For example, the 
orbit and rotation are both assumed to be perfect clocks and, as a result, data time series are often 
” tuned” to the astronomical time scale (Martinson et al., 1982; Shackleton et al., 1990; Hilgren, 
1991). However, as we shall see in the next two sections, the rotational variations in particular 
do not keep very good time. Also, a number of significant observations remain without adequate 
explanation. One of the most perplexing of these is the mid-Pleistocene climatic switch in dominant 
oscillation frequency. The 41 kyr oscillation was dominant throughout the early Pleistocene, and 
the 100 kyr oscillation has been dominant since the Brunhes-Matuyama magnetic polarity reversal 
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(or thereabouts), whereas the standard July 65° N insolation curve shows no discernible change in 
spectral composition over the entire interval. 

The intent of this article is to point out the need for a critical examination of this standard scenario, 
with the objective of improving the geodynamic component of the model enough to provide better 
radiative forcing time series to the paleoclimate community, use paleoclimate proxy data records to 
calibrate solid earth responses to applied loads and torques, and explore potential feed-back loops . 

Two of the topics of study are perturbations to the simple precession model presented in Equation 
(4). The first topic concerns time variations in the precession rate occasioned by glacial mass 
transport and resulting changes in the difference between polar and equatorial moments of inertia 
(C-A). These fractional variations can likely reach 1% in magnitude and are fully competitive with 
changes in orbital eccentricity in terms of their effect on instantaneous precession rate. The second 
topic is the effect of differential precession in the deep interior, which is governed by the (poorly 
known) characteristics of the inertial and dissipative coupling torques which attempt to keep the 
rotation axes of the mantle and core aligned. 

The third topic is the potentially significant geodynamic feed-back loop associated with deglacial 
mass transport, polar motion and ensuing perturbations to radiative equilibrium temperature pat- 
terns. The asymmetrical disposition of major ice sheets and ocean basins relative to the rotation axis 
implies that during growth or decay of these ice sheets the geographic location of the principal axis 
of greatest inertia will shift. If the hydrospheric and cryospheric changes are appreciably more rapid 
than can be compensated by asthenospheric flow, the rotation axis will shift, possibly by as much 
as " 1°. The primary climatic consequence is a shifting of the geographic distribution of continental 
and oceanic regions (which differ by a factor of 60 in heat capacity). The effect of a 1° change in 
the geographic position of the pole is different (but not necessarily less significant) than a 1° change 
in obliquity. 

The fourth and final topic is a new modelling strategy in which climate proxy data records are 
used to estimate linear combinations of insolation pattern Fourier-Legendre coefficients which best 
duplicate the observed variations. This approach allows the data variance to be partitioned into 
global versus regional effects, and can distinguish between responses to annual average insolation 
and those due to seasonal cycle fluctuations. 
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PRECESSIONAL DYNAMICS WITH VARIABLE RATE FACTOR 


Statement of Problem 

All but the most recent reconstructions of the radiative forcing input to paleoclimate models have 
assumed that both the orbital and rotational dynamics could be readily and accurately reconstructed 
from their present configurations, via the simple analyses mentioned in the introduction. These 
expectations seem well founded in the case of orbital evolution, though the possibility of chaotic 
dynamics in the inner solar system (Laskar, 1990; Laskar et ah, 1992) does seem to preclude confident 
extrapolation beyond 10 7 years. However, there are a number of processes, working in different 
locations and at different rates, which all serve to compound the difficulty of accurately computing 
the spin precessional evolution. 

An important aspect of this problem is the synergism between geodynamics and paleoclimate. 
Better geophysical models are needed in order to accurately reconstruct the radiative input to 
climate models. Some of the paleoclimate proxy records contain information relevant to solid Earth 
processes on time scales which are difficult to constrain otherwise. 

On the longest time scales of interest (10 7 -10 9 years) the limiting uncertainty is variability in the 
tidal transfer of angular momentum from the rotation of the Earth to the orbit of the Moon. At 
present, these tidal torques are increasing the length of the day by 22.5 10 6 sec/year and increasing 
the size of the lunar orbit by 3.88 cm/year (Cazenave and Daillet, 1981; Christodoulidis et al., 
1988). Berger et al. (1989) have made a useful first step towards including this effect in climatic 
time series. They computed the change in the major precession and obliquity frequencies due to 
lunar tidal evolution, assuming that the present rate of tidal energy dissipation is representative of 
the past 500 Myr. However, the present rates are considerably higher than the long term average 
(Hansen, 1982), largely due to a near resonance between sloshing modes of ocean basins and the 
diurnal and semidiurnal tidal periods (Platzman et al., 1981), and apparently compounded by a 
contribution from shallow seas (Wunsch, 1986; Dickman and Preisig, 1986). Sedimentary records 
which constrain lunar orbital evolution show some promise of resolving this problem (Olsen, 1986; 
Williams, 1989a, b; Herbert and D’Hondt, 1990), but the situation is definitely more complex than 
is suggested by Berger et al., (1989). 

Another parameter which can vary, on rather shorter time scales and in an equally irregular 
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fashion, is the gravitational oblateness of the Earth (C-A)/C. Thomson (1990) has recently made 
three important contributions to the understanding of this source of variability. First, he pointed 
out that mass redistribution associated with major glaciations and compensating subsidence and 
crustal deformations (Le Treut and Ghil, 1983; Wu and Peltier, 1984) can cause fractional changes 
in oblateness of order 10 ^-10 Second, he showed that high resolution spectral analyses of 
several climatic time series appear to indicate fluctuations of the luni-solar precession rate of this 
magnitude, and with a dominant period near 100 kyr. Finally, Thomson pointed out that the 
best fit to the paleoclimate proxy data was obtained using a mean lunisolar precession rate 0.6 
arcsec/year less than the present observed value. He notes that the resulting value would correspond 
rather closely with that expected for a hydrostatic flattening (Nakiboglu, 1982). If these important 
results are corroborated, they will demonstrate that important feed-back loops exist in the orbital- 
rotational-climatic interactions system, further up-stream” in the presumed causal chain than has 
been previously recognized. 

Approach 

The research design for this segment of the proposed investigation will address several issues. The 
primary focus will be an attempt to resolve two related questions. What are the precessional and 
paleoclimatic consequences of small (0.001-0.01) variations in oblateness (C-A)/C, over time scales 
of 10 -10 years? Can such variations can be confidently inferred from paleoclimatic proxy data? 
These questions represent a forward modeling problem and a coupled inverse problem, respectively. 

The first question is the easier to answer. Equation (4) describes the variations in orientation of 
the spin axis, as viewed in an inertial reference frame. However, since we are primarily interested in 
the orientation of the spin axis relative to orbit normal (obliquity) and the apsidal line (longitude 
of perihelion), the analysis will be made easier if we first transform to a coordinate system fixed in 
the orbit plane. If the rotation matrix is denoted by A, the transformed equation takes the form 
(Ward, 1974; Bills, 1990a) 


ds/ dt = {dA/dt A 1 }s + a(n ■ s)(s x n) 

where 


(9) 


{dA/dt A~ 1 } = B dl/dt + C dW/dt 


( 10 ) 



(11) 


B = 


0 0 0 

0 0 1 

0-10 


0 cos (/) — sin(/) 

C = I —cos (/) 0 0 I (12) 

sin(/) 0 0 

Now define two complex quantities: P = sin(I) e‘ n , which represents components of the orbit 
normal on the invariable plane, and S = sin(c) e’* , which represents components of the spin vector 
on the orbital plane. In this new notation, Equation (9) can be rewritten in the form 


dS/dt + ia S = ib dP/dt (13) 

where 

a = acos(e ) + cos(I) dQ/dt (14) 

b = cos(e)e~ in (15) 

The complete solution to a non homogeneous linear differential equation consists of both ” free” 
and ” forced” modes of oscillation. The free modes, in this case, correspond to spin precession with 
the orbit plane fixed, and the forced modes correspond to motions of the spin axis as it attempts to 
precess about the instantaneous orbit normal, while the orbit normal itself is precessing. The forced 
modes make first order contributions to both c and zv, whereas the free modes are only second order 
for e but are first order for zj. As a result, in the standard linear series solutions (Berger, 1978; 
Ward, 1974), the obliquity terms include only the forced response to changes in (1,0), whereas the 
nodal longitude terms include both forced modes from (1,0), and free modes with variable (e,u/). 

The result of variable oblateness will have exactly the same qualitative effect on spin precession as 
does a change in orbital eccentricity. Both effect the free modes only. It will thus be rather difficult 
to confidently distinguish oblateness variations from unmodeled eccentricity variations. However, 
that does not effect the forward modeling aspect of the problem. To the extent that oblateness 
variations occur, their effect should be included in the astronomical forcing to climate models. 

An iterative approach to the problem seems promising. Orbital variations are uneffected by 
oblateness and need not directly concern us. As a first step, standard rotational variations (including 
eccentricity, but neglecting oblateness variations) will be used to generate radiative input to an 
energy balance climate model (North et al., 1981; Short et al., 1991), with a coupled ice-sheet model 
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(DeBlonde and Peltier, 1990). Time variations in oblateness can be simply estimated from the 
surface load and internal compensation. The resulting oblateness history is then fed back into the 
rotational calculation and the entire process is repeated. The inner-most loops of this algorithm are 
somewhat similar to the model of Peltier (1982). However, he did not allow the mass loading to 
influence the radiative forcing. 

Even fairly modest changes in oblateness are rotationally significant. For comparison, the eccen- 
tricity perturbation influences precession rate via the factor (a/b) 3 , where a and b are semimajor 
and semiminor axes, respectively. This amounts to (1-e 2 ) -3 / 2 , which differs from unity by only 
0.0054 for a near-maximum value of e = 0.06, 
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DIFFERENTIAL PRECESSION: INERTIAL AND DISSIPATIVE COUPLING OF 

THE MANTLE AND CORE 


Statement of Problem 

The hydrostatic figure of a planet represents a compromise between gravitation, which attempts 
to attain spherical symmetry, and rotation, which prefers cylindrical symmetry. Due to its higher 
mean density, the core of the Earth is more nearly spherical than the mantle. The direct lum-solar 
precessional torques on the core will thus be inadequate to make it precess at the same rate as the 
mantle. In fact, the core oblateness is only about 3/4 that required for coprecession with the mantle 
(Smith and Dahlen, 1981). However, it is clearly the case that the core and mantle precess at very 
nearly the same rate (Stacey, 1973). A variety of different physical mechanisms contribute to the 
torques which achieve this coupling, but a purely phenomenological partitioning is useful. The net 
torque can be described as a sum of inertial torques, which are parallel to ( X m * Xc), and dissipative 
torques, which are parallel to ( Xm - Xc)- Here, *c and X m are the rotation vectors of the core and 
mantle, respectively. The two types of torques have qualitatively different results: inertial torques 
cause the core and mantle axes to precess at fixed angular separations and on opposite side of their 
combined angular momentum vector, whereas the effect of dissipative torques is to reduce the angle 

between the axes. 

On short time scales it is appropriate to consider the core to be an inviscid fluid constrained 
to move within the ellipsoidal region bounded by the rigid mantle (Poincare, 1910; Toomre, 1966; 
Voorhies, 1991). The inertial coupling provided by this mechanism is effective whenever the ellipticity 
of the container exceeds the ratio of the precessional to rotational rates. If the mantle were actually 
rigid, or even elastic (Merriam, 1988; Smylie et al„ 1990), this would be an extremely effective type 
of coupling. However, on sufficiently long time scales, the mantle will deform viscously and can 
accommodate the motions of the core fluid (Wu, 1990). The inertial coupling torque exerted by the 
core on the mantle will have the form 

Ti = ki[ Xm x Xc] (16) 

A fundamentally different type of coupling is provided by electromagnetic or viscous torques 
(Rochester, 1962; Sasao et ah, 1977; Kubo, 1979). The dissipative coupling torque exerted by the 
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core on the mantle will have the form 


Td — Xc] 


(17) 


This type of coupling is likely to be most important on longer time scales. In each case, the mantle 
exerts an equal and opposite torque on the core. The response of the coupled core-mantle system 
to orbital forcing is given by (Goldreich and Peale, 1970; Ward and DeCampli, 1979; Bills, 1990b). 


dsm/dt — a m (n s m )(s m x n) — 0 m {s m — s c ) — 7 m (s m x s„) 

( 18 ) 

ds c /dt = a c (n ■ s c )(s e x n) + /3 c (s m - s c ) + y c (s m x s c ) 
where o m is similar to o above, except that only mantle moments A m and C m are included, and 


Pm 

7m 

w T here v is the mean rotation rate. 


kdj Crn V 

ki/C m v 2 


(19) 


Approach 


The research design for this segment of the proposed investigation would consist of several parts. 
The objectives and methods are very similar to those described in the previous section. In this 
case, however, the intent would be to determine the precessional and paleoclimatic consequences of 
non-rigd core-mantle coupling, and to explore the possibility that useful constraints on the coupling 
parameters can be obtained from paleoclimatic proxy data. 

A number of estimates already exist for the strengths of inertial and dissipative coupling torques 
(Toomre, 1966, Roberts, 1972; Rochester, 1976; Loper, 1975; Stix, 1982). By most accounts, the 
inertial torque is ' 5 10 20 N m, and the various viscous and magnetic dissipative torques are 102-104 
times weaker. However, the inertial torque estimates are simply based on the premise that the core 
must cop recess with the mantle. 

Solutions to the coupled precession problem can be found in a form analogous to Equation (11) 
for the rigid precession problem. The principle difference in the present situation is the increased 
richness of the free and forced oscillation spectra. There are modes in which the core and mantle 
precess together, and other modes which reflect differential precession. It is clear that the most 
climatically relevant behavior is the precessional motion of the mantle. Thus the chief interest, from 
that perspective, will be to explore the behavior of the mantle precession modes over plausible range 
of parameter values. Within the geophysically prescribed range of parameter space, are there any 
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climatically significant perturbations to obliquity or longitude of perihelion occasioned by partially 
decoupling the core from the mantle? 

It is very clear that the mantle and core exhibit differential motions on nutational time scales 
(Mathews et al., 1991; Mathews and Shapiro, 1992). However, one of the difficulties in constructing 
a differential precession model which is truly useful for paleoclimatic studies is that the precise 
geodetic techniques, which are able to constrain nutation amplitudes within a few milliarcseconds, 
still have short enough time spans that most of the differential precession modes of interest are still 
indistinguishable from rigid rotation. 
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DEGLACIAL POLAR MOTION AND ABRUPT CLIMATE CHANGE 


Statement of Problem 

Recent models of the surficial mass transport associated with the last deglaciation (Tushingham 
and Peltier, 1991; Nakada and Lambeck, 1988) suggest that the magnitude of the flow and its 
departure from axial symmetry were both great enough that, if it were not closely balanced by 
internal flow and deformation, it would cause a shift in the body-fixed location of the axis of greatest 
inertia, by an amount of order lo. If the pole were to move by a significant fraction of a degree, the 
resulting displacement of the geographic pattern of land and water relative to the incident radiation 
pattern would cause a climatic perturbation which could either augment or retard the progress of 
deglaciation, depending on the spatio-temporal pattern of the perturbation. The potential thus 
exists for a previously unexplored feed-back loop in the climate system. 

During the peak of the last deglaciation, there was a brief but significant return to full glacial 
conditions. Ibis Younger Dryas climatic event is perhaps best documented in the North Atlantic 
(Ruddiman and McIntyre, 1981; Broecker et al., 1988), but appears to have been global in scale 
(Currey, 1990; Engstrom et ah, 1990; Gasse et ah, 1991). Though the broad scale timing of the 
deglaciation is consistent with astronomical forcing, the Younger Dryas perturbation was too brief to 
be a linear response to the classical orbital or rotational variations. However, neither the magnitude 
nor duration is a priori inconsistent with a response to deglacial polar motion. The objective of this 
portion of the proposed investigation would be to examine the geodynamic and climatic consequences 
of deglacial mass flow. 

Approach 

The usual geodynamic modeling approach to deglacial polar motion is to use historical observations 
of the rate and direction of polar motion during this century, in conjunction with estimates of the 
surficial mass flow, to derive constraints on deep Earth structure models. However, these same 
models, once they have been calibrated, can be made to deliver estimates of the rates and directions 
of polar motion which occurred during the deglaciation. The magnitude of deglacial polar motion 
depends rather sensitively on the rate and spatial pattern of surficial mass transport and on the 
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internal structure of the Earth, as reflected in the spectrum of relaxation times for surficial loads 
and body forces (Sabadini et ah, 1982; Wu and Peltier, 1984; Spada et ah, 1992; Ivins et ah, 1992). 

A class of climate models which is very well suited to the proposed investigations incorporates 
spatial and temporal variations in a number of components of the global energy balance. The basic 
equation which governs the energy balance is (North et ah, 1981) 

A+ BT + C^-V {DVT) = QaF (20) 

at 

The first two terms on the left parameterize outgoing radiation, the third term represents storage 
of heat and the fourth term represents divergence of the heat flux. The right hand side is the amount 
of incident radiation that gets absorbed. 

Values for the outgoing radiation parameters (A = 210 W m -2 , B = 2.1 W m K 1 ; Short 
et ah, 1984) can be estimated from satellite data. The heat capacity C is large over water and 
small over land (C„,/B = 4.6 years, C, = C„, /60; North and Coakley, 1979). Though much of the 
lateral transport of heat is accomplished by winds and ocean currents, these effects can be modeled 
by a diffusive transport, with D/B = 0.310 (North, 1975; Wyant et ah, 1988). For the incident 
radiation, Q = S/4 = 342 W m -2 (Schatten and Orosz, 1990), F is the spatial and temporal pattern 
of projected area, and a (0.75, Stephens et ah, 1981) is the co-albedo. 

Solutions to this equation with specified forcing can be readily obtained in the spectral domain, 
with spatial patterns represented by spherical harmonic series and temporal patterns represented by 
Fourier series. If A, B and D are assumed to be constants (no spatial or temporal variations), the 
base state solution to Equation (20) is simply 

Tnmp = W~ l [Q G a nmp F n0p ] (21) 

where G (whose 9 indices have been suppressed for typographical clarity) is a coupling constant for 
products of spherical harmonics (Rotenberg et ah, 1959; Dill, 1991) and 

W = [B + n(n+ 1) D+ipG Cnmp] (22) 

In the same notation, the temperature perturbation 6T due to a small change in radiative forcing 
<$F (associated with c, e, or w) is 

8Tnmp = W- 1 [Q G a n mp 6F n0p ] (23) 

Similarly, the temperature perturbation due to a change <5C in the distribution of land and water 
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( 24 ) 


relative to the rotation axis (associated with deglacial polar motion) is 

fiT nmp — W 1 [—ipG SC nmp T nmp ] 

Specific questions which should be addressed include: 

• What is the likely history (magnitude and direction) of polar motion induced during deglacia- 
tion? 

• How does the climatic impact of that polar motion compare in magnitude and spatio-temporal 
pattern with the impact of a 1° change in obliquity, or a 1% change in orbital eccentricity? 

• Is the timing and magnitude of the climatic impact such that it would significantly perturb 
the deglaciation? 

• Is this feed-back loop a viable contributor to the Younger Dryas climatic event? 
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FOURIER, LEGENDRE, AND MILANKOVITCH: ESTIMATING DIAGNOSTIC 
SEASONS AND LATITUDES FROM CLIMATE PROXY RECORDS 

Statement of Problem 

As was mentioned in the introduction, a weakness of previous paleoclimate models is the com* 
mon practice of using an overly simplistic representation of the insolation forcing (values for July 
at 65° N are a particular favorite) when attempting to reconcile computed variations in orbital and 
rotational parameters with observed climate proxy records (Hays et ah, 1976). An alternative ap- 
proach, which shares many of the same problems is to use a linear combination of the orbital and 
rotational parameters as input to the models (Imbrie and Imbrie, 1980). A problem in using either 
a linear combination of orbital and rotational parameters, or a localized spatio-temporal measure 
of insolation, as a proxy input to climate models is the difficulty in choosing which single value to 
use (Broecker, 1966; Broecker and Van Donk, 1970). Furthermore, the two alternative formats for 
making the choice (time and location versus orbital and rotational parameter combination) are not 
equivalent. It is clear that choice of a specific latitude and time of year to represent the insola- 
tion pattern implies a unique (though non-linear) combination of orbital and rotational parameters. 
However, the converse is not necessarily true. For some combinations of (c,e,ru) there is no corre- 
sponding location and time of year, at which the insolation pattern will be dominated by the selected 
combination of parameters. 

In comparing calculated insolation pattern variations with an observed time series of some climate 
proxy (marine sediment oxygen isotope variations, for example), it is useful to have a simple physical 
model in mind. The observed isotopic anomalies are, to first order, due to two effects; global ice 
volume fluctuations and local to regional scale temperature variations (Emiliani and Shackleton, 
1974; Kahn et al., 1981). These effects differ in two important ways. The ice volume effect is 
global in scale and its first time derivative should be proportional to insolation driven temperature 
changes. In contrast, the direct isotopic temperature effect is local to regional in scale, and is directly 
proportional to insolation driven temperature changes. 
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Approach 


A modeling strategy which addresses several of these points involves estimating the linear com- 
bination of Fourier-Legendre coefficients of the insolation pattern which best reproduces the data 
records. This approach provides the opportunity to partition observed variations into global and 
regional effects. It also provides information on which aspects of the insolation pattern variations 
are most diagnostic of the observed proxy variations. If July insolation at 65° N is truly signifi- 
cant in a particular data record, the selected amplitudes for the Fourier-Legendre coefficients should 
clearly reflect that fact. Alternatively, if the amplitude of the seasonal cycle in the polar regions or 
tropics, or the annual mean equator to pole insolation contrast, or any other linear functional of the 
insolation pattern is most diagnostic, the analysis will indicate that fact. 

As an example of a situation in which this approach could be used, Park and Maasch (1992) 
have recently compared the climate proxy record provided by <$ 18 0 records from two long, high 
resolution cores (ODP 677 from the eastern equatorial Pacific (Shackleton et ah, 1990) and DSDP 
607 from the mid-latitude North Atlantic (Ruddiman et al., 1989; Raymo et al., 1989) with a 
new estimate of insolation at 65° N (Berger and Loutre, 1991). The three data records analyzed 
(benthic and planktonic records at ODP 677, benthic only at DSDP 607) have much in common, 
but there are also interesting, and possibly significant, differences. It is clear that variations which 
are common to several locations are more likely diagnostic of global climatic variations. Differences 
between locations, or differences between benthic and planktonic records at a single location provide a 
different and complementary view of the climatic response to insolation forcing. Figure (5) illustrates 
a 2 10 6 year history of variations in the amplitudes of the Fourier-Legendre coefficients and Figure 
(6) compares the DSDP 607 and ODP 677 data with a standard insolation curve. 

Several specific questions should be addressed. The first category of questions relate to calibration 
and characterization of technique. 

• Do the Fourier-Legendre coefficient time series form an orthogonal basis? 

• What output time series would correspond to a spatio-temporal delta function (summer solstice 
at 65° N, for example)? 

• What is the resolution (A/i,AM) versus data series length? 

• How does the resolution depend on other factors besides the data string length? 
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• Are two short data strings from significantly different locations better than a single data string 
of equivalent aggregate length? 

The next level of question pertains to the physical meaning or significance of a particular linear 
combination of Fourier-Legendre coefficients. Any such linear combination is equivalent to a spatio- 
temporal pattern of some sort. The result obtained with a finite resolution data structure will 
represent the convolution of the actual physical response with the data resolution kernel. How 
do errors in data chronologies map into errors in resolved spatio-temporal patterns of diagnostic 
insolation? 

The final level of question relates to understanding the implications and significance of results 
obtained by this algorithm from actual data series. 

• Using a variety of oceanic and continental paleoclimatic data series, what are the resolved 
patterns of climatically diagnostic insolation? 

• What physical mechanisms are implicated? 

• Are significantly different patterns obtained from data strings of different age? 

» Is there a single pattern which duplicates the observed increase in amplitude of the 100 kyr 
signal after the Brunhes-Matuyama magnetic reversal? 
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Figure 2 

Orbital inclination spectrum and time series. 

Both are based on the secular variation model of Laskar (1988). 
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Figure 3 

Rotational obliquity spectrum and time series. 

Time series is based on rotational theory of Kinoshita (1977) and orbital model of Laskar (1988). 
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Figure 4 

Latitudinal and Seasonal Insolation Pattern. 
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Figure 5a 

Insolation Pattern Fourier-Legendre Coefficient Time Series. 
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Figure 5b 

Insolation Pattern Fourier-Legendre Coefficient Time Series. 
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Figure 6 

Paleoclimate Proxy Records and Insolation Time Series. 

Curves are (top to bottom): DSDP 607 (benthic), ODP 677 (planktonic), ODP 677 (benthic), 65° N insolation 
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Introduction 

The insolation parameters of the Earth depend on its orbital parameters and of the pre- 
cession and obliquity. Until 1988, the usually adopted solution for paleoclimate computation 
consisted in (Bretagnon, 1974) for the orbital elements of the Earth, which was completed by 
(Berger, 1976) for the computation of the precession and obliquity of the Earth. In 1988, I 
issued a solution for the orbital elements of the Earth, which was obtained in a new manner, 
gathering huge analytical computations and numerical integration (Laskar, 19SS). In this 
solution, which will be denoted La88, the precession and obliquity quantities necessary for 
paleoclimate computations were integrated at the same time, which insure good consistancy 
of the solutions. Unfortunately, due to various factors, this latter solution for the precession 
and obliquity was not widely distributed (Berger, Loutre, Laskar, 19SS). On the other side, 
the orbital part of the solution La88 for the Earth, was used in (Berger and Loutre, 1991) 
to derive an other solution for precession and obliquity, aimed to climate computations. I 
also issued a new solution (La90) which presents some slight improvements with respect to 
the previous one (Laskar, 1990). As previously, this solution contains orbital, precessional 
and obliquity variables. In order to make it widely available, it was distributed during this 
meeting on magnetic support and can also be obtain directly by request to the author at 
laskarOf riap51 .bitnet. In the present talk, I will discuss the main features of this new 
solution. 
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The orbital solution La90 

The orbital solution is obtained by the numerical integration of an extended averaged system, 
which represents the mean evolution of the orbits of the planets. All the S main planets of the 
solar system are taken into account, as well as lunar and relativistic main perturbations. The 
use of numerical integration for the computation of the solution of the secular system is one of 
the reasons for the good quality of this solution, which can be checked by comparison with the 
available ephemeris over a short time scale (Laskar, 1985, 1986, 1988). In (Laskar, 1988), the 
solution LaSS was represented in quasi-periodic form over 10 Myr, but these representations 
are slowly convergent, which prevents a good accuracy for the solution. Later on, I understood 
fully the reason of this slow convergence, which is due to the presence of multiple resonances 
m the secular system of the inner solar system. Due to these resonances, the motion of the 
solar system is chaotic, and not quasi-periodic, as was demonstrated by the computation of 
its Lyapunov exponents which reachs 1/(5 Myr) (Laskar, 1989). This implies that it is not 
possible to give any precise solution for the motion of the Earth over more than about 100 
Myr, and most probably, ephemeris can only be given for about 10 Myr with good precision. 
Several integrations of the secular system of the solar system were made over 200 Myr and 
400 Myr. The origin of the chaotic behaviour was identified, and is due to the presence of 
secular resonances in the inner solar system (Laskar, 1990). With a new numerical method, 
it was possible to show that the chaotic zones were the solar system belongs is large in the 
directions of the proper modes related to the inner planets. This is an indication that these 
results are stable against small changes of initial conditions or model. 

One of the main consequence of interest for paleoclimate studies, is the fact that the main 
requencies of the orbital motion of the planets can no longer be considered as constants, but 
are slowly evolving with time. The measured shift in frequency amount about 0.2 arcsec/year 
over 200 Myr for ^3 and g 4 while for g 5 it is of only 0.00002 arcsec/year. It should be stressed 
that, as the motion is chaotic, the computed solution cannot be considered as close to the 
real solution over more than a few 10 Myr, but it is reasonable to assume that the drift 
in frequency is of the same order of magnitude over the 200 Myr. One should nevertheless 
mention that over this time span, the change of the precession main frequency, due to tidal 
eftects in the Earth- Moon system, are probably more important (Berger et ai, 1989) 

Since, direct numerical integrations over 3Myr backward and forward were issued by (Quinn 
et ai, 1991), including also solutions for precession and obliquity (QTD6). The orbital 
solutions have been compared with La90. The two solutions present very small discrepancies 
over 3Myr, and the existence of the secular resonances in the inner solar system is confirmed 
( askar, Quinn, Tremaine, 1991). The very close agreement of the two orbital solutions over 

.* ' yr glVeS the Confirmatlon that th e Earth parameters are now very well known over this 
time span and insure that the orbital solution La90 can be used with confidence over 10 Myr 

or pa eochmate use. The precession and obliquity solutions present some small discrepancies 

^ h . ar ' P robabl y due to the presence of the tidal effect in the Earth-Moon system in the 
Q 1D6 solution. 
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Precession and Obliquity in the La88 and La90 solutions 

The precession quantities are completely determined by the two motions of the equatorial 
and ecliptic pole. The motion of the ecliptic is given by the secular theory La90 (Laskar 
1990); the precession quantities are integrated at the same time, using the equations of the 
theory of the rigid Earth of Kinoshita (Kinoshita 1975, 1977, Laskar, 1986). The equations 
for the general precession in longitude p.4, and for the obliquity of the date e are then 


with : 


and : 


dpA 
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where p = sin(r‘/2) sin(fi), q = sin(i/2) cos(fl), (i is the inclination of the Earth with respect 
to a fixed ecliptic, and Q the longitude of the node). R(e) is the secular term due to the 
direct lunisolar perturbations. The quantities Mo, Mi, M2, M3, So, and S2 depend only on the 
orbital elements of the Moon and the Sun. The pincipal moments of inertia of the Earth are 
denoted by A, B , and C, and the angular velocity of the Earth is 10 . The masses of the Sun, 
the Earth, and the Moon are denoted by m©, and ; the sideral motion of the Sun and 
of the Moon by n© and riM ; and the mean motion of the node of the Moon by uq. The other 
terms present in equation (23) represent the effects of the secular variation of the ecliptic, 
caused by the secular planetary perturbations. The numerical values of Mo, Mi, M2, M3 are 
given in (Kinoshita, 1977): 

Mo = 496303.3 x 10~ 6 
Mi = -20.7 x 10 -6 
M 2 = -0.1 x 10 -6 
M 3 = 3020.2 x 10 -6 
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and from (Laskar, 1986), 


So 


i(l -e 2 ) -3 / 2 - 0.422 x 10 -6 


The following numerical values, were also used (see Laskar, 1986 for complete references) : 


w = 474 659 981 .597 57 "/yr 

n M = 17325593.4318 "/yr 

n Q = -69 679 .193 6222 "/yr 

a M = 384 747980.645 m 

k = 0.01720209895 


m©/(mE + itim ) = 328 900.5 
m © l m E — 332 946.0 


The quantity p g is the geodesic precession due to the general relativity, 


p g = 0.019188'Vyr 

The value of the dynamical ellipticity E D = (2 C -A- B)/2C = 0.0032S005 is obtained by 
adjustment at the origine J2000 to the values of the speed of precession and obliquity give 
by the IAU (Grenoble, 1976) : 

p = 50.290966"/yr 
e 0 = 23°26'21"448 

For t = 0 , we have p A = 0 ,£ = e 0 , i = = 0 ,and thus : 


dpA 

dt 


= R{to ) - 2 p t=0 - p g cot e 0 

t=o 


These formula for precession gives a solution for precession and obliquity in agreement with 
the requirements of high accurate ephemeris (Laskar, 1986) and are thus best suited for 
paleoclimate computations. 


Descripion of the files of the solution La90 

The different files which are distributed are 
ear0m5.dat orbital elements 0 to -5 Myr 

pre0m5.dat obliquity and precession 0 to -5 Myr 

ear5ml0.dat orbital elements -5 to -10 Myr 

pre5ml0.dat obliquity and precession -5 to -10 Myr 

All the elements are referred to the ecliptic and equinox J2000. The starting date is J2000. 
All angles are in radians 

The ear0m5.dat and ear5ml0.dat files contains T, k, h,q,p 
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T = time from J2000 in lOOOyr 
k = ecos(ru) 
h = e sin(c7) 
q = sin(i/2) cos(fl) 
p = sin(z/2) sin(fi) 

where e is the eccentricity of the Earth, i denotes the inclination of the Earth, w the longitude 
of perihelion, and Q the longitude of the node of the Earth with respect to the fixed ecliptic 

and equinox J2000. 

The pre0m5.dat and pre5ml0.dat files contains T,e,p A where 

T is the time from J2000 in lOOOyr, £ is the obliquity of the date (mean equator of date with 
respect to mean ecliptic of date), and p A is the general precession in longitude. 

All quantities generally used for climate studies are derived easily from these fundamental 
quatities. The eccentricity of the Earth is obtained by 

e = %A 2 + V 

And the longitude of perihelion from moving equinox is 

U* = ZD + PA 

The climatic precession e cos(u>*) is thus equal to the real part of z*exp(ip A ), where 2 = k + ?h. 

Conclusions 

The present solution La90 for orbital and precession elements for the Earth over can be used 
for as an ephemeris paleclimate computations over 10 Myr. On this time span, the fact that 
the motion of the solar system is chaotic is not perceptible on this level of precision. This 
new solution present some improvements from my previous solution LaSS (Laskar 19S8). The 
orbital solution is in very good agreement with the recent numerical integrations of (Quinn et 
al., 1991) . Over the geological time scales exceeding 100 Myr, there is no hope to obtain such 
an ephemeris, due to the chaotic behaviour of the solar system, but one can assume that for a 
few 100 Myr, the slow diffusion of the astronomical frequencies observed during the 200 Myr 
integrations remains of about the same order. This should not be granted for billion years 
computations and more computations on the diffusion process in the solar system are clearly 
needed. More, due to the changes in frequencies, and presence of secular resonances in the 
orbital forcing, more extended studies of the resulting influence on the rotational evolution of 
the Earth should be done. Studies on the long term evolution of the Earth-Moon system are 
also needed in order to improve our knowledge of the rotational evolution of the Earth. The 
geological records, assuming that a good accuracy in the determination of the fundamental 
frequencies could be achieved, would then be the only possible observations for tracking the 
long term past evolution of the solar system for time span longer than 100 Myr. 
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1 Overview 

A theory for the long-term rotational motion of the quasi-rigid Earth has been constructed 
by numerical integration. The theory spans 72000 years centered about 1968 A.D., and 
provides accurate rotational and positional data for the Earth in the recent past and the 
near future. Details are provided in [6] . The physical model is termed dynamically consistent 
because developments for the active forces and torques are truncated based solely on their 
magnitudes regardless of their origin. The model includes all appropriate forces and torques 
due to the geopotential and tidal effects as well as lunisolar and planetary contributions. 
The elastic and inelastic deformations due to tidal action were too small to affect the mass 
properties of the Earth at the truncation level of the model. However, long-term dissipative 
effects of the tidal forces and torques were not negligible. These considerations gave the 
model its quasi-rigid characterization. The numerical output provided both rotational and 
orbital-element data. The data have been fitted throughout the 72000-year range using 
Chebyshev polynomial series. These series sire quite portable and are available upon request. 

The project was based on the desire to provide researchers with a data base of ac- 
curate orbital and rotational parameters which are needed as the astronomical input to 
(paleo)climatology theories. Current theories [2] use only a rigid-Earth model and are 
analytical in nature. In addition, they are not completely consistent in their physical mod- 
els. To keep pace with the accuracy of current observational data, the numerical theory 
maintains an accuracy on the order of 0.01 mas. The physical model includes: 

1. rigid-body torques produced by the Moon, Sun, Venus, Jupiter and Saturn 

2. tidal torques due to the Moon and Sun 

3. effects due to Earth rotation/lunar orbit coupling 

4. point-mass effects arising from the orbital motion of a 10-body planetary system 
(counting the Sun and Moon with Pluto excluded) 

5. a 4x4 geopotential field 
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The numerical integration was conducted using a 10 th -order Adams-Moulton predictor 
with an ll th -order Adams-Bashforth corrector. This combination produced a local 
truncation error E « /i 12 . The stepsize h was set to 1/50 day. The selection of h was 
based on error comparisons with output generated by extended precision (32 place) 
integrations. 

The major results are displayed graphically in [6]. The data give an average precession 
rate of -50.45' ; /year. The present adopted LAU value is — 50.29"/year, and that of 
the best analytical theory ([2]) is -50.4l"/year. The output for obliquity (relative to 
the invariable plane of the solar system) and the precession index 1 esinu; are shown 
on the following page. 


2 Future Work 

The data base is not nearly long enough to be of great use to paleoclimatologists. Part 
of what is needed as an accurate theory that extends backward approximately 1 million 
years. (The data would thus span at least two of the 400,000-year Milankovitch cycles) . 
A hybrid numerical integration procedure would be required — one that effectively 
suppresses the growth of round-off error. Such a procedure has been developed by 
Panovsky and Richardson [4] and has been used successfully in a multi-million-year 
integration of the planetary system [5]. Data spanning a much longer time frame 
would allow the retrieval of an accurate (geologically-recent) history for the lunar 
retreat rate. 

To produce a data base spanning a much greater time frame, the quasi-rigid model 
for the Earth in the present theory would have to be discarded in favor of some sort 
of radially-dependent viscoelastic model which includes core-mantle coupling (gravi- 
tational and pressure) and mass-property variations. The new model should include 
changes in continentality from plate motions as well as the long-term effects of the 
oceans and glaciation. It may be possible to devise a model that is accurate over the 
long-term by adjusting parameters so that the (concocted) tidal dissipation provides 
lunar retreat rates consistent with the various measurements that have been reported 
in the recent literature. The extremes of these measurements now span 650 million 
years (see [l] and [7]). 

To push a numerically-integrated theory backwards through tens of millions of years 
would require the development and use of a semi-numerical theory applied to a much- 
improved physical model. Such a theory is based on the numerical integration of 
the equations of motion after they have been analytically averaged at least through 
second order. A similar process has been successfully implemented by Laskar [3] in his 
recent investigations of the long-term motion of the eight-body planetary system. The 
averaging process would remove all high-frequency effects and leave only the secular, 
long-period, and dominant non-linear effects. The numerical integration could then 
proceed with stepsizes on the order of weeks or months instead of hours. 

lr The parameters c and Cj are respectively, the orbital eccentricity of the Earth and the Earth’s argument 
of perihelion measured from the moving Vernal equinox. 
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Variations of the Milankovitch Frequencies in 

Time. 


Loutre M.F. and Berger A. 


1 Pre-Cenozoic Times 


The sensitivity of the amplitudes and frequencies in the development of the Earth’s 
orbital and rotational elements involved in the astronomical theory of paleoclimates 
(eccentricity, obliquity and climatic precession), to the Earth-Moon distance and con- 
sequently to the length of the day and to the dynamical ellipticity of the Earth has 
been discussed for the last billions of years (Berger et al., 1989a, b,c, Berger and Loutre, 
1991). 

The shortening of the Earth-Moon distance and of the length of the day, as well 
as the lengthening of the dynamical ellipticity of the Earth back in time induce a 
shortening of the fundamental astronomical periods for precession (the 19-kyr and 23- 
kyr quasi-periods becoming respectively 12.6 and 14.3-kyr at 2 10 9 yr BP) and obliquity 
(the 41-kyr and 54-kyr quasi-periods becoming respectively 19.6 and 22.1-kyr at 2 10 9 
yr BP) (Figure 1). At the same time, the amplitudes of the different terms in the 
development of the obliquity are undergoing a relative enlargement of about 50% at 
2 10 9 yr BP but the independent term is increasing very weakly (less than 0.1%). In 
other words, the value of the obliquity, which lies within a range of 21. 7 to 24. 9 
over the Quaternary was restricted to a range of 22. °5 to 24.° 1 at 2 10 yr BP. On 
the other hand, the amplitudes in the development of the climatic precession do not 
change. Moreover, these changes in the frequencies and amplitudes for both obliquity 
and climatic precession are larger for longer period terms. Finally, the periods in the 
eccentricity development are not influenced by the variation of the lunar distance. 


But the motion of the solar system, especially of the inner planets, was shown to 
be chaotic (Laskar, 1990). It means that it is impossible to compute the exact motion 
of the planets over more than about 100 Myr, and the fundamental frequencies of the 
system are not fixed quantities, but are slowly varying with time. As long as we consider 
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the most important terms, the maximum deviation from the present-day value of the 
19-kyr precessional period due to the chaotic motion of the solar system only does not 
reach more than a few tens of years around 80 Myr BP (Berger et al., 1991). Therefore 
the shortening of the obliquity and climatic precession periods is mostly driven by the 
change in the lunar distance and the consequent variations in the dynamical ellipticity 
of the Earth’s angular speed. 

At first sight, the deviation in the period for the eccentricity can be neglected, as the 
chaotic behaviour of the solar system implies a relative change of the main periods (404, 
95 and 123 kyr-periods respectively) by less than 0.2%, 1.4% and 1.9% respectively, 
this maximum changes being achieved around 80 Myr BP. This implies, in particular, 
that the eccentricity periods used for Quaternary climate studies may be considered 
more or less constant for pre-Quaternay times and equal to their Quaternary values. 


2 The Quaternary Period 


The sensitivity of the frequencies of these astronomical elements to the dynamical 
ellipticity of the Earth has also been investigated for the Quaternary period (Dehant 
et al., 1990). According to the model used, the modification of the distribution of the 
masses on and inside the Earth during full glacial conditions has a weak influence on 
the moments of inertia of the Earth and consequently on the astronomical periods: 
more precisely the dynamical ellipticity of the Earth remains more or less constant for 
the isostatic equilibrium case (i.e. a subsidence of about 3/10 of the height of the ice 
sheet). For larger subsidences, the variation of the dynamical ellipticity of the Earth is 
positive and for smaller ones, it becomes negative, reaching 1.5% at the maximum. For 
a 1 % increase of the dynamical ellipticity of the Earth, the precessional quasi-periods 
(19 and 23 kyr) become 18.9 and 22.8 kyr while the quasi-periods of the obliquity (41 
and 54 kyr) become 40.4 and 52.9 kyr. 


3 Future research 

• The determination of the Earth-Moon distance and/or the lunar recession rate 
for pre-Cenozoic times must be improved by taking into account the repartition 
and the displacement of the oceanic basins and of the continents in order to 
obtain a better time scale. 

• A full model of the Earth interior would give more accurate values for the dy- 
namical ellipticity of the Earth and the Earth’s angular speed of rotation taking 
into account the repartition of the masses on and inside the Earth. 

• The model used to compute the variation of the dynamical ellipticity of the Earth 
due to the waxing and waning of ice sheets could be improved by considering the 
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— geographical location and height of the ice sheets 

— subsidence and rebound of the continents not covered by ice 

— subsidence and rebound of the marine lithospere 

• A model accounting for the transient response to the ice sheet loading, instead 
of the snapshot reconstruction as used in the present study where we considered 
only full (maximum) glacial conditions and plain interglacials, would allow us to 
simulate the time evolution of the global effects of the ice sheets formation and 
melting, taking into account lagging subsidence and rebound of the lithosphere. 

• Improvements in the interpretation of past proxy data and refinement of the time 
scale would allow to provide an independent determination of the astronomical 
frequencies found in pre-Cenozoic time series, allowing to validate and calibrate 
the astronomical models. 
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Figure 1: Estimated values of the periods of the orbited parameters (top) involved 
in the astronomical theory of paleoclimates considering the effect of the variations 
of the Earth-Moon distance (bottom) and of the Earth’s figure and rotation (Berger 
and Loutre, 1991). It must be noted that the discontinuity in the rate of change of 
the astronomical periods reflects the artificial change in the value of the Earth-Moon 
recession rate taken to be 10 -9 m s -1 for the last 590 Myr , and 0.43 10 -9 m s -1 prior 
to 590 Myr BP (Walker and Zahnle, 1986; Berger and Loutre, 1991). 
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The long, well-exposed and often continuous sequences of loess found throughout the world 
are generally thought to provide an excellent opportunity for studying long-term, large-scale 
environmental change during the last few million years. In recent years, the most fruitful loess 
studies have been those involving the deposits of the loess in China. One of the most intriguing 
results of that work has been the discovery of an apparent correlation between variations in the 
magnetic susceptibility of the loess sequence and the oxygen isotope record of the deep sea. This 
correlation implies that magnetic susceptibility variations are being driven by astronomical 
parameters. However, the basic data have been interpreted in various ways by different authors, 
most of whom assumed that the magnetic minerals in the loess have not been affected by post- 
depositional processes. Using a chemical extraction procedure that allows us to separate the 
contribution of secondary pedogenic magnetic minerals from primary inherited magnetic minerals, 
we have found that the magnetic susceptibility of the Chinese paleosols is largely due to a 
pedogenic component which is present to a lesser degree in the loess. We have also found that the 
smaller inherited component of the magnetic susceptibility is about the same in the paleosols and 
the loess. These results demonstrate the need for additional study of the processes that create 
magnetic susceptibility variations in order to interpret properly the role of astronomical forcing in 
producing these variations. 

The Chinese loess plateau stretches from 35°N to 40°N and from 100°E to 1 15°E and covers an 
area of 500,000 sq. km. The loess deposits are typically 150 m thick, and they appear to represent 
continuous deposition of wind-blown, silt-sized material during the past 2.4 million years. The 
source of this material is believed to be glacial outwash in the regions to the west and north of the 
plateau (Kukla and An, 1989). More importantly, the loess sequence contains many interbedded 
paleosols which attest to the existence of significant and cyclic climatic fluctuations. The most 
recent, comprehensive description of the units of the loess sequence is that of Kukla and An 
(1989), who recognized six stratigraphic units. From youngest to oldest, these are the Holocene 


Black Loam Formation, the Malan Formation, the Upper Lishi Formation, the Lower Lishi 
Formation, the Wucheng Formation, and the Pliocene Red Clay. 

On the basis of paleomagnetic studies (Heller and Liu, 1982; Kukla, 1987), the contact 
between the Wucheng Formation and the Red Clay layer has been dated at 2.4 million years, and 
the Brunhes/Matuyama, Olduvai event and Jaramillo event have each been identified in the 
sequence. 

One of the primary parameters that has been used in the study of the loess/paleosol sequence 
has been magnetic susceptibility. Because the magnetic susceptibility of the loess is low while that 
of the paleosols is high, this parameter is considered an effective proxy for the quantitative study of 
the climatic fluctuations recorded by the loess/paleosol sequence (Heller and Liu, 1984; 1986). 

The first comprehensive study of magnetic susceptibility variations in the loess/paleosol 
sequence was that of Heller and Liu ( 1984) who pointed out that there appeared to be a strong 
correlation between the magnetic susceptibility record and the oxygen isotope record of deep-sea 
cores from the equatorial Pacific Ocean. The relationship was further explored by Kukla et al. 
(1988) who published detailed magnetic susceptibility records from the loess/paleosol sections at 
Xifeng and Luochuan. These authors presented data to support their belief that the time required 
for the deposition of a particular loess unit was directly proportional to the product of the thickness 
of the unit and its magnetic susceptibility. They used this idea to construct a time scale that was 
independent of the oxygen isotope curve. On this time scale, the variations in magnetic 
susceptibility corresponded very closely to the variations in the oxygen isotope record from the 
deep sea, implying an interdependence among the rate of influx of loess, the volume of land-based 
ice, and the global climate. Additional evidence for astronomical forcing of the magnetic 
susceptibility record was provided by Wang et al. (1990). 

A key component in the model used by Kukla et al. (1988; 1990) to account for the magnetic 
susceptibility variations was the assumption that the source of the magnetic susceptibility signal 
was a constant "rain" of ultrafine magnetic grains, carried into the upper atmosphere from volcanic 
eruptions and other unspecified processes. Kukla et al. further assumed that after these grains had 
been incorporated into the loess sequence during deposition, they remained inert and unaltered by 
post-depositional processes. The loess, on the other hand, was assumed to be essentially non- 
magnetic, and the modulation of the magnetic susceptibility signal was interpreted as a measure of 
the extent to which the magnetic "rain" had been diluted by loess. Thus, during glacial times, 
when the climate was cold and dry, the barren outwash plains could be easily eroded by aeolian 
processes, the rate of loess deposition would be at a maximum, and the magnetic susceptibility 
signal would be at a minimum. During interglacial times, when the climate was warm and humid, 
vegetation and soil moisture would tend to stabilize the outwash plains, loess deposition would be 
a minimum, and the magnetic susceptibility would be a maximum. 
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The model of Kukla et al. (1988; 1990) differs from that of Heller and Liu (1984) who 
suggested that the magnetic susceptibility values in the paleosols reflected a concentrating of the 
magnetic minerals by decalcification and soil compaction. Both models discounted any post- 
depositional alteration of the magnetic carriers. This fundamental assumption has been questioned 
by Zhou et al. (1990), Maher and Thompson (1991), and Zheng et al. (1991) who showed that 
there were significant differences between the rock magnetic properties of the magnetic minerals in 
the loess units and those in the paleosol units. These differences implied that there were 
differences in both the magnetic mineralogy and the grain size of the magnetic minerals in the two 
units. Maher and Thompson (1991) also raised questions about the methods that Kukla et al. used 
to demonstrate that the rate of accumulation of magnetic minerals had been constant. Zhou et al., 
Zheng et al., and Maher and Thompson all concluded that pedogenic processes had probably been 
important in the development of the magnetic susceptibility record of the paleosols. 

We have obtained direct evidence that the magnetic susceptibility signal of both the loess units 
and the paleosols is due primarily to magnetic minerals formed by pedogenic processes. This 
conclusion is based on the studies of samples from ten loess/paleosol pairs from the classic section 
in Luochuan. The samples were provided to us by George Kukla of the Lamont- Doherty 
Geological Observatory, and they encompass the entire loess/paleosol sequence. Their 
designations, stratigraphic positions and approximate ages are shown in Table 1. 

For each sample, we measured a variety of rock magnetic properties both before and after 
extraction with citrate-bicarbonate-dithionite (CBD). In this procedure, samples are subjected to 
sodium dithionite and bicarbonate, a strong buffered reductant, in the presence of sodium citrate, a 
chelating agent (Singer and Janitzky, 1986). The procedure was developed by Mehra and Jackson 
(1960) as a means of removing iron oxides from clay samples being prepared for X-ray diffraction 
analysis. The procedure was subsequently adopted by soil scientists as part of the standard 
chemical technique for characterizing the iron components of a soil. With that technique, extraction 
procedures involving pyrophosphate, oxalate and CBD are used to determine the amount of iron in 
organic, amorphous and crystalline phases, respectively. In recent years, we have used the CBD 
extraction technique in our studies of magnetic susceptibility enhancement in soil chronosequences 
in California (Singer and Fine, 1989; Fine et al., 1989; Singer et al., 1992). That work has shown 
that that CBD extraction is particularly effective in removing pedogenic magnetic grains (primarily 
maghemite) and that it leaves untouched essentially all of the magnetic grains that were inherited 
from the soil parent material (primarily magnetite and hematite). This selectivity has recently been 
confirmed by Mossbauer spectrometry (Singer et al., 1991). 

For untreated samples from the loess plateau, our rock magnetic measurements are fully 
consistent with those reported by Maher and Thompson (1991) and by Zhou et al. (1990). For 
example, the magnetic susceptibilities of the paleosols are as much as twenty times larger than 


51 


those of the corresponding loess samples (Table 1). Differences between loess and paleosol 
sample pairs are also noted in the frequency dependence of the magnetic susceptibility (which is a 
measure of the concentration of ultrafine grained, supeiparamagnetic particles), in the S-ratio 
(which is a parameter related to hematite concentration), in the ratio of saturation isothermal 
remanent magnetization to anhysteretic remanent magnetization (which is a measure of the relative 
abundance of single domain grains) and in the ratio of magnetic susceptibility to anhysteretic 
remanent susceptibility (which is related to mean magnetic grain size). 

After CBD treatment, both the loess and the paleosol samples lose a significant percentage of 
their magnetic susceptibility (Table 1). These losses average 65% for the loess samples and 90% 
for the paleosol samples. Because the magnetic susceptibilities of the untreated paleosol samples 
are five to ten times greater than that of the untreated loess samples, the absolute decreases in 
magnetic susceptibility are much greater in the paleosols than in the loess units (Figure 1). 
Furthermore, after CBD treatment, the magnetic susceptibilities of the loess samples and the 
paleosol samples are about the same, regardless of the age of the samples (Figure 1). Several other 
rock magnetic properties also show decreases after CBD treatment with largest changes again 
occurring in samples from the paleosols (Figure 2). For a few rock magnetic properties, the values 
from the paleosol and loess samples move in opposite directions after CBD treatment. 

Furthermore, we have found a close relationship between magnetic susceptibility and dithionite- 
extractable iron (Figure 3), providing additional evidence of the importance of pedogenesis in 
determining the magnetic susceptibility of paleosols and loess. 

Based on our work on the soil chronosequences in California, we interpret the CBD soluble 
fraction in the loess and paleosol samples as the pedogenic fraction, and the CBD insoluble fraction 
as the inherited fraction. This indicates that a significant portion of the magnetic susceptibility 
signal of both the loess samples and the paleosol samples is pedogenic in origin. The fact that 
pedogenesis is important in producing the magnetic susceptibility signal in the paleosols was 
suggested by Zhou et al. (1990), Maher and Thompson (1991) and Zheng etal. (1991). 

However, none of these groups postulated that pedogenesis could account for almost all of the 
magnetic susceptibility signal in the paleosols, and none of them proposed that pedogenesis would 
be important in the loess units as well. This latter observation gives us an entirely new perspective 
on the paleosol/loess sequences. In the conventional view, paleosol units are considered to have 
resulted from very different processes than those that produce the loess units. From our results, it 
seems that the same pedogenic processes might have been operating during times of loess 
deposition and paleosol formation but these processes were more intense during the former than 
during the latter. 

Our results also show that other earlier inferences about the nature of the magnetic 
susceptibility signal were probably also wrong. For example, Maher and Thompson suggested 
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that the pedogenic component of the magnetic susceptibility was probably carried by magnetite 
while our data strongly support the conclusion that maghemite is the primary magnetic mineral. 
More importantly, our observations provide no support for the concept of an inert, ultrafine 
magnetic "rain" diluted to varying degrees by non-magnetic windblown silt, as proposed by Kukla 
et al. (1988; 1990). In fact, if the nearly constant residual magnetic susceptibility that we observe 
in both the paleosol and loess samples after CBD treatment is an exogenous magnetic component, 
it implies that the loess was accumulating at a same rate during glacial and interglacial stages and 
that the differences between paleosols and loess are due entirely to the degree of pedogenesis. 

At present time, we are not prepared to argue the merits of this or any other explanation of our 
results. What we will argue is that we have shown that there is a clear need for a better 
understanding of the nature and origin of the magnetic susceptibility signal in the Chinese 
loess/paleosol sequence. This need is more than just a minor problem, of interest to a small group 
of rock magnetists. As noted above, the loess/paleosol sequences in general, and the Chinese 
sequences in particular, are considered the best recorders of terrestrial climate change during the 
last 2.4 million years. Almost exclusively, this change is being studied using magnetic 
susceptibility as a proxy indicator of paleoclimate. In fact, using assumptions about the magnetic 
susceptibility signal that our research has now shown to be incorrect, other workers have already 
developed an elaborate model for climate changes in Asia and the western Pacific. The model 
attributes these changes to astronomically-driven fluctuations in the summer monsoon that are 
modulated by uplift of the Tibetan plateau (Kukla, 1987; An et al., 1991). 

While certain aspects of this model may ultimately prove to be correct, the model itself cannot 
be validated until its underlying assumptions are based on the proper paleoclimate interpretation of 
the magnetic susceptibility record. Our work has shown that this interpretation must address the 
pedogenic nature of the magnetic susceptibility signal. This requirement also applies to 
loess/paleosol sequences elsewhere that are also being interpreted as records of terrestrial climate 
change. 
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Table 1. 


Paleosol and Loess Samples from Luochuan, China. 
(Samples provided by George Kukla.) 


Position 

Depth 

Aget 

Magnetic susceptibility 
Pre-CBD Post-CBD 


m 

ka 

SI x 10-8 

kg" ' 

S 1 

10.0 

1 28 

227.5 

22.7 

L2-LLI 

12.5 

1 74 

1 1.4 

8.1 

S3 

24.0 

3 28 

223.3 

19.7 

L4 

26.0 

357 

57.2 

18.6 

S 5 

36.5 

6 14 

283.3 

1 1.4 

L6 

42.0 

652 

59.7 

13.7 

S 7 

52.0 

726 

81.1 

15.3 

L8 

52.5 

7 

58.5 

15.0 

S 8 

54.5 

737 

129.5 

15.0 

L9 

57.5 

834 

26.9 

15.3 

L14 

76.0 

? 

53.7 

1 1 .2 

L 1 5 

78.5 

1 172 

23.0 

14.4 

WS1WLI 

86.0 

13 16 

132.9 

13.9 

WS1SS1 

86.5 

? 

86.2 

14.3 

WS1SS1 

87.5 

? 

84.6 

14. 1 

WS2SS2 

99.0 

1566 

102.5 

15.0 

WL2LL2 

101. 0 

1695 

54.9 

12.8 

WS3LL1 

107.5 

1939 

45.0 

13.0 

WS3LL1 

109.0 

1939 

49.9 

9.4 

WL4LL3 

134.0 

2342 

37.1 

1 1 .2 

RsSSl 

136.0 

>2342 

168.7 

14.3 


t Age from Kukla (1987) Table 5. 
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Figure 1 Effect of CBD treatment on magnetic susceptibility of some paleosol and loess units 
from the Lishi Formation, Luochuan, China. Squares are paleosol units; circles are loess units 
Arrows indicate change upon CBD treatment. The post-CBD values for both paleosols and loess 
units are about the same. 
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Frequency dependent susceptibility 
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Figure 2. Effect of CBD treatment on the frequency dependence of the magnetic susceptibility of 
some paleosol and loess units from the Lishi Formation, Luochuan, China. Squares are paleosol 
units; circles are loess units. Arrows indicate change upon CBD treatment. The post-CBD values 
for both paleosols and loess units are about the same. 
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Figure 3. Relationship between dithionite extractable iron and magnetic susceptibility of some 
paleosol and loess units. 
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1 Climatic CoolinP related to G eomagnetic Reversals! 

There have been speculations on the relationship between climatic cooling and 
polarity reversals of the earth's magnetic field during the Pleistocene (Kawai et al. 1975, 
Rampino 1981, Krishnamurthy et al 1986, and Jacobs 1984 for a review). Two of the 
common criticisms on this relationship have been the reality of these short duration 
geomagnetic events and the accuracy of their dates. Champion et al. (1988) have reviewed 
recent progress in this area. They identified a total of 10 short-duration polarity events in 
the last 1 Ma and 6 of these events have been found in volcanic rocks, which also have K- 
Ar dates. The nomenclature and the estimated ages for these events are shown on the top 
of Fig. la. Events with age dated from volcanic rocks are shown as solid bars, those from 
sedimentary rocks as stippled bars. Also, typeface and size of the names represent the 
degree of confidence in the exact age or existence of that polarity event. Following 
Rampino (1981), the eccentricity of the earth's orbit over the past 1.2 Ma is also calculated 
(Berger 1977) and plotted with the stacked oxygen isotope data (with the SPECMAP time 

scale, Imrie et al 1984) at the bottom of Fig. la. 

An inspection of Fig. la shows that during the last 600 ka, where the 100 ka cycle 

(due to eccentricity) was the dominant signal in the oxygen isotope data, magnetic polarity 
events seem to occur near times of maximum eccentricity and rapid glaciation, thus giving 
some support to a possible climate-magnetic reversal connection. For t > 600 ka, the 
correlation is not as good - although there is an apparent correlation between times of 
maximum eccentricity and the Kamikatsura event, the end of the Jaramillo and the Cobb 
Mountain event. If the ages of these older events are revised according to Shackleton et al. 
(1990), then all of these older events, except the Kamikatsura, appear to occur near 
eccentricity maxima (Fig. lb). 

Anyhow, Champion et al. (1988) found that the mean of the polarity interval 
lengths (in Fig. la) to be close to the 100 ka main orbital eccentricity period of the earth, 
they therefore suggested that linkage between geomagnetic, paleoclimatic and possible 
underlying earth orbital parameters should be further studied. 
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Fig.la Correlation of the polarity events with the eccentricity of the Earth’s orbit, and 
80 18 (ice volume) variation over the past 1.2 Ma. 



Fig.lb Same as Fig.la, except ages greater than 600 ka are modified according to 
Shackleton eL al (1990) 
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2. Mechanisms that relate climatic cooling & Geomagnetic Rgygrsals 

Supposing that the speculated relationship between climatic cooling and geomagnetic 
reversals actually exist, two mechanisms that assume climatic cooling causes short period 
magnetic reversals will be investigated. It should be noted that this is NOT an attempt to 
explain the occurrence of ALL the magnetic polarity events in terms of glacial advances, for 
it is obvious that magnetic reversals has occurred throughout the history of the earth - even 
during times when there is no glaciation. Therefore, the mechanisms that we investigate 
operates only within the Pleistocene and it should be clear that other mechanisms, with 
different time scales, may be operative at the same time and they may be responsible for the 
reversals outside this epoch. 

2.n Core-Mantle Boundary Topography 

A possible mechanism results from the variation of topographic interaction across the 
core-mantle boundary (CMB). For example, the formation of large bumps or depressions 
at the CMB will set up Taylor columns in the core (Hide, 1969), disrupt the flow which 
drives the geodynamo and possibly change the geomagnetic field. The critical height for 
topographic coupling to be effective has been estimated to be around 1 km . 

Gubbins and Richards (1986) have investigated the effect of thermal and subduction- 
induced topography and have concluded that both have the right magnitude to cause 
coupling. However, the thermal time scale is of the order of 100 Ma whereas changes in 
topography due to subduction is about 10 Ma. Since polarity reversals can occur several 
times in 1 Ma, there must be other phenomena which cause the short time scale variations 
in topographic coupling during the Pleistocene. Now, the period of a glacial cycle is about 
100 ka and for a 3 km thick ice sheet, the depression at the earth's surface is about 1 km. 
Thus, depending on how rapid this deformation attenuates with depth, glacial induced 
topography at the CMB, if suitably located, may be able to enhance (or diminish) the 
coupling due to thermal and subduction-induced topography such that the total coupling can 
exceed (or fall below) the threshold value required to disrupt core flow. 

The topography of the CMB due to glacial loads at the earth's surface has been 
calculated for 3 different earth models LI, L2 & L3 which has lower mantle viscosities of 
1021, 3xl()21 \ q 22 p a . s respectively. The ice model consists of 3 centers of 

glaciation corresponding to the Laurentide, FennoScandian and Antarctic ice masses which 
has total mass equivalent to a sea level drop of 100 meters and whose glacial history 
consists of 30 cycles of glaciation, the last of which is plotted in Fig. 2a. The result of this 
calculation is described in Wu (1990) which shows that the maximum topography occurs 
underneath the Laurentide ice center and the time variation of this maxima is shown in 
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Fig.2c. In view of its smallness in height and its location relative to other topographic highs 
as deduced from seismic tomography, glacial induced CMB topography is not expected to 
be able to significantly modulate total coupling and modify core flow today. 

Transfer of Rotation al Energy to the core 

An alternate mechanism is the transfer of rotational energy into the fluid core: the 
redistribution of water masses during glaciation and deglaciation will cause the moment of 
inertia of the earth to change - for example, when water is taken from the ocean basins and 
is locked in the ice sheets near the pole, the moment of inertia will decrease. By the 
conservation of angular momentum, the decrease in the moment of inertia must be 
accompanied by an increase in the angular velocity of the mantle. If the increase in velocity 
at the bottom of the mantle exceeds the 0.03 cm/sec flow velocity in the fluid core 
(estimated from the westward drift of the non-dipole component), then the transfer of 
angular momentum into the core would disrupt the convective heat engine inside and 
possibly result in a change in magnetic polarity (Doake 1977, Muller & Morris 1986). 

The change in velocity at the bottom of the mantle due to the simple 
glaciation/deglaciation model described earlier is plotted in Fig.2b. The solid line, the 
dotted line and the dashed line correspond to earth models LI, L2 and L3 respectively. The 
stripped area in the middle indicates that the velocity at the bottom of the mantle is below 
the preexisting flow velocity in the fluid core. From Fig.2b it can be seen that model LI 
(and possibly L3 near glacial minimum) can produce velocity changes comparable to the 
existing flow in the core. 

Given that a 100m change in sea levels can produce velocity changes at the bottom of 
the mantle to be comparable to the existing flow in the core, the next question concerns 
how energy is transferred from the mantle to the core. The mantle and the core can interact 
by electromagnetic coupling, viscous coupling and topographic coupling. If the roughness 
of the CMB exceeds the thickness of the viscous hydromagnetic boundary layer, then 
topographic coupling becomes the dominant mechanism . Although the topography induced 
by glacial loads has been shown to be too small to cause coupling, the topography due to 
density loads in the mantle exceeds the critical value and possibly provides the coupling 
mechanism to transfer rotational energy into the core to cause geomagnetic reversals. 

3. Conclusion 

In conclusion, the variation in CMB topography induced by the surface glacial loads 
has the correct time scale but not the amplitude, unless the maxima are suitably located, to 
significantly modify the total coupling mechanism and disrupt the flow pattern in the core 
resulting in geomagnetic field reversals. The transfer of rotational energy from the mantle to 
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the core via topographic coupling is a more likely candidate, provided that the change in sea 
level is of the order of 100 meters and the coupling is via thermal and subduction induced 

topography. 

4. Questions and Future Work 

More precise dates and the establishment of the reality of some polarity events (eg. 
the Emperor & Delta events) are required to establish/disprove the controversial relationship 
between climatic cooling and geomagnetic reversals. In the compilation of Champion et al 
(1988), some of the events have rather large uncertainties in their age, and, as shown in 
Fig.l, the correlation of some of the older events with eccentricity maxima depends 
critically on their age. 

Opponents to a climate- magnetic reversal relationship often question whether polarity 
events observed in different parts of the world but with approximately the same age can be 
correlated and whether they correspond to the same global events. Another question is 
whether there are only 8 polarity events in the Brunhes, for, if more events are discovered, 
then there may be no correlation between polarity events and eccentricity maxima. Clearly, 
cores in different parts of the earth, with more complete record and dates are needed to 

answer these questions. 

Even when the relationship between climate and magnetic reversals is established, the 
next question is which is the cause and which is the effect? If magnetic field is the cause, 
then the question is why they do not always cause climatic cooling (eg. during the 
Cretaceous)? If climatic cooling is the cause, then, are there other mechanisms? 

A better seismic tomographic map of the CMB would confirm whether glacial 
induced topography can modulate total coupling and modify core flow. Glacial induced 
topography is rejected because, in seismic tomographic maps of the CMB (Morelli & 
Dziewonski 1987, Creager & Jordan 1986), no topographic hills/troughs seem to exist 
underneath Laurentia. These seismic tomographic maps are, however, not consistent with 
each other and therefore have some degrees of uncertainty in them. 

More work is needed to understand the details of the core-mantle coupling 
mechanisms and how energy and angular momentum can be transferred from the mantle to 

the core. 

Finally, even after this is achieved, questions still remain as to how this energy and 
momentum are going to disturb the flow field and the convective engine of the core? how a 
polarity reversal comes about ? and why the reversals are often so brief? 
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Fig.2 (a) The last cycle of the saw-tooth load history, (b) The time history of the change in 
velocity at the bottom of the mantle due to increased spin rate of the mantle, (c) The time 
evolution of the CMB topography under the Laurentide ice center for the 3 earth models. 
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Late Quaternary Time Series of Arabian Sea Productivity: 
Global and Regional Signals 


S.C. Clemens W.L. Prell D.W. Murray 


Brown University 
Providence, RI 


Modern annual floral and faunal production in the northwest Arabian Sea derives primarily from 
upwelling induced by strong southwest monsoon winds during June, July, and August. Indian Ocean 
summer monsoon winds are, in turn, driven by differential heating between the Asian continent 
and the Indian Ocean to the south. This differential heating produces a strong pressure gradient 
resulting in southwest monsoon winds and both coastal and divergent upwelling off the Arabian 
Peninsula (Figure 1). Over geologic time scales (10 4 to 10 6 years), monsoon wind strength is sensitive 
to changes in boundary conditions which influence this pressure gradient. Important boundary 
conditions include the seasonal distribution of solar radiation, global ice volume, Indian Ocean sea 
surface temperature, and the elevation and albedo of the Asian continent. To the extent that these 
factors influence monsoon wind strength, they also influence upwelling and productivity. In addition, 
however, productivity associated with upwelling can be decoupled from the strength of the summer 
monsoon winds via oceanic mechanisms which serve to inhibit or enhance the nutrient supply in the 
intermediate waters of the Indian Ocean, the source for upwelled waters in the Arabian Sea (Prell, 

1990). 

To differentiate productivity associated with wind-induced upwelling from that associated with 
other components of the system such as nutrient sequestering in glacial-age deep waters (i.e. Boyle, 
1988) we employ a strategy which monitors independent components (Figure 2) of the oceanic and 
atmospheric subsystems. Using sediment records from the Owen Ridge, northwest Arabian Sea, we 
monitor the strength of upwelling and productivity using two independent indicators, % G. bulloides 
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(Prell, 1984, 1990) and opal accumulation (Murray, 1990). We monitor the strength of southwest 
monsoon winds by measuring the grain-size of lithogenic dust particles blown into the Arabian Sea 
from the surrounding deserts of the Somali and Arabian Peninsulas (Clemens and Prell, 1990a). 

1 he planktonic foraminifer G. bulloides is typically a subpolar species found in the Southern ocean 
between the Subtropical Convergence and the Antarctic Convergence. However, G.bulloides is also 
found in high abundance in tropical upwelling areas such as the northern Arabian Sea of the coast of 
Arabia (Hutson and Prell, 1980; Cullen and Prell, 1984) and Cariaco Trench, located off Venezuela 
(Overpeck et al., 1989). G. bulloides abundance in the Arabian Sea is negatively correlated to 
summer sea surface temperature (Prell and Curry, 1981) which, in turn, is negatively correlated 
with wind-induced upwelling during the summer monsoon (Prell and Streeter, 1982). Similarly, opal 
accumulation is positively correlated with nutrient distribution associated with upwelling but less 
influenced by regional sea surface temperatures (Murray, in preparation). 

The Owen Ridge lies beneath the axis of the strong summer monsoon winds. These winds can 
transport lithogenic dust particles such as those found on the Owen Ridge (mean diameter of 14.4 
fim) for thousands of kilometers given mean velocities of 15 m/s and values for the coefficient of 
turbulent exchange found in typical cyclonic storms and frontal systems (Tsoar and Pye, 1987). 
The largest particles (40 to 50 /xm) can also be transported up to 1000 km in more extreme dust 
storm events. These estimates of transport distances are consistent with studies which identify the 
Somali and Arabian Peninsulas as primary source areas for dust found in Arabian Sea sediments 
(Middleton, 1986; Sirocko and Samthein, 1989). Over geologic time scales, increases in the strength 
of monsoon winds result in the transport of larger lithogenic particles to the Arabian Sea, thus 
increasing the median grain size of the lithogenic component (Clemens and Prell, 1990a). 

Spectral analyses of these records allow us to examine concentrations of variance held in common 
between these independent abiotic and biotic records of wind strength and productivity. Our results 
unambiguously demonstrate that all three independent records are linearly related and in phase 
with one another over the Earth’s orbital precession cycles (23 kyr cycles; Figure 3). We interpret 
these relationships as indicating that: (1) to first order, all three indicators are linked by a common 
response to monsoon wind strength and upwelling, and (2) precessional insolation is the primary 
external forcing mechanism for the late Quaternary monsoon and the associated upwelling induced 
productivity. However, variance associated with precession accounts for only " 25% of the total 
variance in any given record. A large portion of the total variance in the biotic records is concentrated 
at the 41 kyr period associated with the obliquity of the Earth’s orbit (Figure 4; Prell, 1990). The 
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grain-size record, on the other hand, contains large amounts of variance at the 29, 35, and 54-kyr 
periods which represent heterodyne periods of the primary orbital periods (100, 41, and 23 kyr 
periods; Clemens and Prell, 1990b). It is the dissimilarities between these records that contains 
the information possibly allowing differentiation of productivity associated with monsoon-induced 
upwelling from that associated with other mechanisms. 

Our current hypothesis is that the variability associated with the 41 kyr power in the G. bulloides 
and opal accumulation records derive from nutrient availability in the intermediate waters which 
are upwelled via monsoon winds. This hypothesis is testable by comparison with Cd records of 
intermediate and deep waters of the Atlantic and Indian Ocean (e.g. Boyle, 1988). 

The 35 and 54 kyr heterodyne periods in the grain-size record are positively correlated with 
large amplitude insolation events at 30°south in the Indian Ocean, the latitude of the subtropical 
high pressure cells from which the Indian Ocean monsoon winds initiate. The 29 kyr variability 
is linearly related to sea surface temperature records from the southern subtropical Indian Ocean. 
Both associations can be explained via the relationship between latent heat (a function of ocean- 
atmosphere temperature gradients) and monsoon strength as follows (Clemens and Prell, 1990b). 
Latent heat collected overt the southern subtropic Indian Ocean is transported across the equator 
and released in the mid-troposphere about the Tibetan Plateau. I his increases the strength of 
the monsoon low, resulting in increased wind strength and transport of larger lit hie particles to 
the Arabian Sea. This hypothesis is currently being tested by development of late Quaternary sea 
surface temperature records from ~ 30° south. Concentrations of variance at the 35, and 54 kyr 
periods in these records would support the latent heat link between the insolation record and the 
grain size record of monsoon strength. 

Confirmation of the hypotheses described above will eventually enable us to quantitatively par- 
tition variance within these records into that attributed to monsoon-induced upwelling (regional) 
and oceanic mechanisms of nutrient distribution (global). Identification of a global signal in regional 
upwelling productivity would then provide a framework for comparison of productivity records from 
upwelling regions throughout the world ocean. 
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Location Map 
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winds (arrows) flow from high pressure to low pressure inducing coastal and divergent upwelhng 
(shaded) off the Arabian Peninsula. Upwelling productivity is recorded m the geological record of 
fossil planktonic foraminifera (carbonate) and radiolaiia (opal) preserved in Owen g 
sediments. Both southwest and northwest summer winds transport terrigenous du « jo tiifr Owen 
Ridge from deserts of Somali and Arabia. The grain size of the lithogemc component vanes as a 
function of the strength of the transporting winds. 
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Figure 2. Lithogenic and biogenic indicators of monsoon wind strength and upwelling-induced 
productivity. Increased % G. bulloides and opal accumulation record changes in productivity 

TnrSS? WUh nutne ^ t , co " t ? nt of * e intermediate waters upweUed during the summer monsoon. 
“ S 1 " 310 Slze °f the hthogenic component varies as a function of the strength of summer 
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Figure 3. Precessional phase wheel summarizing the coherency and phase relationships between 
several monsoon indices and the Earth's orbital precession. All parameters shown are coherent (at 
or above the 80% confidence level) with insolation patterns driven by the precession of the Earth's 
orbit The phasing indicates that all the monsoon records attain maxima ~9 kyrs after maxima in 
precessional insolation but simultaneously with minima in sea surface temperatures (SST) at 13° 
south in the Indian Ocean. This, similar to modem monsoon dynamics, indicates that latent heat 
availability in the southern subtropic Indian Ocean exerts a strong influence on the timing of strong 
monsoons and upweliing in the Arabian Sea over the past 400 kyrs. 
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Frequency Domain 



Hgure 4. Spectra of the biogenic and abiogenic indicators of monsoon wind strength and the 
associated up welling- induced productivity. The 41 kyr (orbital obliquity) variance in the G 
bulloides and opal records is absent from the lithogenic (abiotic) record indicating that variance in 
this frequency band may not be driven by productivity due to monsoon-induced upwelling. This 
variance may be associated with more global oceanic mechanisms which enhance or reduce the 
nutrient content of the intermediate waters of the Indian Ocean. The 29, 35 and 54 kyr spectral 
peaks in the grain size record are linearly related to SST and insolation patterns in the southern 
subtropic Indian Ocean indicating a link between latent heat flux and monsoon strength as is 
observed in modern monsoon dynamics. 
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Annual, Orbital, and Enigmatic Variations in Tropical 
Oceanography Recorded by the Equatorial Atlantic Amplifier 
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Equatorial Atlantic surface waters respond directly to cahuges in zonal and meridional lower 
tropospheric winds forced by annual insolation. This mechanism has its maximum effect along the 
equatorial wave guide centered on 10 -W. The result is to amplify even subtle tropical cl, mate changes 
such that they are recorded by marked amplitude changes in the proxy signals. Model realizations, 
NCAR AGCM and OGCM fo, 0 K, and 126 Ka (January and July), and paleoc.anograph.c proxy 
data show that these winds are also forced by insolation changes at the orbital periods of Precession 


and Obliquity. 

Perhelion in boreal summer produces a strengthened monsoon, e.g. increased meridional and 
decrease zonal wind stress. This reduces oceanic Ekm.n divergence and the, mocline/nutr, dine 
shallowing. The result, in the equatorial Atlantic, is reduced primary productivity and higher 
euphoric zone temperatures; vice versa fo, perihelion in boreal winter. Perihelion is controlled by 
precession. Thus, the dominant period in spectra from a stacked SST record (0-252 Ka BP) at 
the site of the equatorial Atlantic amplifier is 23 Ky (53% of the total variance). This precessional 
period is coherent (k=0.920) and in phase with boreal summer insolation. 

Oscillations of shorter period are present in records from cores sited beneath the amplifier re- 
gion. These occur between 12.5 and 74.5 Ka BP, when eccentricity modulation of precession is 
at a minimum. Within this time interval there are 21 cycles with mean periods of (3.0±0.5)Ky. 
Similar periods have been documented from high latitude regions, e.g. Greenland ice cores from 
Camp Century. The Camp Century signal in this same time interval contains 21 cycles. A subjec- 
tive correlation was made between the Camp Century and the equatorial records; the signals were 
statistically similar, (r=0.722) and (k= 0.960) . 
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Modeling Orbital Changes on Tectonic Time Scales 


Thomas J. Crowley 
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Geologic time series indicate significant 100 ka and 400 ka pre- Pleistocene climate fluctuations, 
prior to the time of such fluctuations in Pleistocene ice sheets. The origin of these fluctuations 
must therefore depend on phenomena other than the ice sheets. In a previous set of experiments 
(Short et al., Quat. Res. 35, 157-173, 1991) we tested the sensitivity of an energy balance model to 
orbital insolation forcing, specifically focusing on the filtering effect of the Earth s geography. We 
found that in equatorial areas, the twice-yearly passage of the sun across the equator interacts with 
the precession index to generate 100 ka and 400 ka power in our modeled time series. The effect is 
proportional to the magnitude of land in equatorial regions. We suggest that such changes may reflect 
monsoon al variations in the real climate system, and that subsequent wind and weathering changes 
may transfer some of this signal to the marine record. A comparison with observed fluctuations of 

Triassic lake levels is quite favorable. 

A number of problems remain to be studied or clarified: 

• the EBM experiments need to be followed up by a limited number of GCM experiments; 

• the sensitivity to secular changes in orbital forcing needs to be examined; 

• the possible modifying role of sedimentary processes on geologic time series warrants consid- 
erably more study; 

• the effect of tectonic changes on Earth’s rotation rate needs to be studied; 

• astronomers need to make explicit which of their predictions are robust and geologists and 
astronomers have to agree on which of the predictions are most testable in the geologic record. 
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Plio-Pleistocene time evolution of the 
100-ky cycle in marine paleoclimate records 

Jeffrey Park and Kirk A. Maasch* 

Department of Geology and Geophysics 
Box 6666, Yale University 
New Haven, CT 06511 


To constrain theories for the dynamical evolution of global ice mass through the 
late Neogene, it is important to determine whether major changes in the record were 
gradual or rapid. Of particular interest is the evolution of the neai 100-ky ice age 
cycle in the middle Pleistocene. We have applied a new technique based on multiple 
taper spectrum analysis which allows us to model the time evolution of quasi-periodic 
signals [Park and Maasch, 1992], This technique uses both phase and amplitude 
information, and enables us to address the question of abrupt versus giadual onset 
of the 100-ky periodicity in the middle Pleistocene. 

The variation of climate proxy variables at a given frequency f a (with associated 
period T 0 = 1 / f a ) can be parameterized by &{ A(t)e~ 2,riJot }, a sinusoid with a slowly- 
varying amplitude .4(f). The Junction ,4(f) is complex-valued, allowing slow variations 
in phase as well as magnitude. We define /1(f) as the envelope of the signal at carriei 
frequency’ f a . If A(t) = A 0 is a constant, the signal is termed ‘periodic’ or ‘a phase- 
coherent sinusoid.’ If A(t) varies, the signal is termed ‘quasi-periodic.’ All methods 
of determining the envelope of a quasi-periodic signal have shortcomings. Bandpass 
filters, such as those used by Ruddiman et al [1989], cannot model discontinuities in 
the envelope of a quasi-periodic signal, and thus, are of little help in discriminating 
between a sudden and a gradual onset lor the 100-ky ice age cycle. Moreover, a 
narrow bandpass in the frequency domain requires a long time-domain filter, so that 
a significant number of data points at the ends of the series must be discarded. 
Similarly, the choice of length for the overlapping data segments in the Ruddiman et 
al [1986b] analysis involves a tradeoff between frequency and time resolution. 

We apply instead a more flexible approach based on multiple-taper spectrum 
analysis [Thomson 1982; 1990; Park et al. 1987, Berger et al 1991]. The multitaper 

‘Current. Address: Institute for Quaternary Studies, .‘120 Boardman Hall, University of Maine, 
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approach allows the analyst to model the envelope function A(t ) using the tools of lin- 
ear inverse theory, solving for the envelope function that fits the time series data while 
optimizing some property of the envelope. Inversion algorithms can be derived that 
find the ‘smallest’ and ‘smoothest’ envelopes that fit spectral ‘data’ from a given time 
series [Park 1990; Park 1992]. Special cases of these the algorithms can be derived 
to model sudden changes in the envelope. This allows us to examine the abruptness 
of the onset of the 100-ky periodicity as well as the evolution of the obliquity and 
precession signals. The estimation procedure relates spectral estimates at / 0 , using a 
set of orthogonal Slepian data tapers, as linear functionals of an unknown envelope 
A(t), so that the envelope is constructed Irom a linear combination of the Slepian 
tapers. Since the Slepian tapers are optimally bandlimited, multitaper envelope es- 
timation can be thought ol as an extremely sharp narrow-band filter valid for the 
entire duration of the series. 

The shortcomings of the algorithm we use are similar to those of other linear 
inveise problems. For instance, the multitaper envelope estimation procedure can 
model envelope discontinuities, but the technique cannot by itself discriminate be- 
tween continuous and discontinuous models for the 100-ky cycle, since envelopes of 
both types can be constructed that fit the spectral ‘data’ exactly. The analyst must 
use other, perhaps subjective, criteria for choosing among models that fit the data. 
This nonuniqueness is a common problem in linear inverse theory; in principle an 
infinite number of envelopes can fit a finite number of spectral data exactly. In 
most cases, the a priori constraint that the signal at a carrier frequency f a have a 
‘slowly- varying’ envelope, except at one or more specified time points, reduces greatly 
the number of acceptable models for the envelope, and motivates the solution of a 
‘smoothest’ envelope that fits the spectral data. 

We investigated the time-evolution of the 100-ky cycle in 6 18 0 data, thought to 
reflect global ice-volume variations. We analyzed 6 18 O data from DSDP Site 607 
and ODP Site 677, from which three long (> 2.6 My) time series have been published 
[Shackleton and Hall, 1989; Ruddiman et al 1986ab; 1989; Raymo et al 1989], We find 
evidence for a coherent 6 18 0 signal for both cores in the eccentricity and obliquity 
frequency bands, consistent with variations in global ice volume as the causative 
factor. However, the nature of the earth system response to orbital insolation cycles 
depends on the time scale adopted in the spectral analysis. If the Ruddiman/Raymo 
time scale for DSDP Site 607 is accepted, the 6 ™ O obliquity cycle has enhanced 
amplitude between 1.0 and 1.5 Ma, relative to the late Pleistocene (t < 1.0 Ma), and 
a nonlinear response of the earth system to orbital obliquity is inferred (Figure 1). 

If the Shackleton et al [1990] time scale for ODP Site 677 is accepted, the amplitude 
match between the 8 18 0 obliquity cycle and the 65°N insolation derived by Berger 


80 



and Loutre [1991] from the recent astronomical solution of Laskar [1988; 1990] is very 
good for times <<2.3 Ma (Figure 2). (The phase match for obliquity is enhanced by 
the fact that the data series were tuned to the astronomical time series, so that the 
observed phase agreement is not surprising.) The veracity of the ODP 677 time scale 
therefore correlates with a linear earth-system response to orbital obliquity. 

Based on our analysis of data from these two deep-sea cores, we do not find 
compelling evidence for an abrupt change in the 100-ky <!> 18 0 signal. Rather, envelope 
inversions in the eccentricity band suggest that the 100-ky 6 18 0 cycle is phase-locked 
with the 124-ky eccentricity cycle some 300-400 ky prior to its late- Pleistocene growth 
in amplitude and phase-lock with the 95-ky eccentricity cycle (Figure 3). An abrupt 
change in the 95-ky envelope near 0.85 Ma is consistent with DSDP 607 data on the 
Ruddiman/Rayrno time scale, but such a transition would occur while leaving the 
124-ky envelope largely unchanged. If the Shackleton et al [1990] time scale for ODP 
677 is accepted, our three 6 18 0 records are consistent with a low-amplitude 100-ky 
cycle between 1.2 and 2.6 Ma, whose local period of oscillation alternates between the 
two major eccentricity periods at 95 and 124 ky. The phase of these 100-ky oscillations 
prior to 1 .2 Ma is related to the phase of the astronomical eccentricity cycles where 
the cycles have significant amplitude. The match between the precessional envelopes 
of Sites 607 and 677 is poor, when both are expressed in terms of the Shackleton et al 
[1990] time scale. Climate simulation studies suggest that cyclic salinity changes in 
equatorial surface waters are a plausible contributor to the ODP 677 <5 18 0 data in the 
precessional band, and could explain this discrepancy. However, our t ime-rescaling of 
data from DSDP 607, using visual isotope-stage matching, may not be precise enough 
to transform the short-period precession cycles with sufficient accuracy. 

Further study of more paleoclimate records will be necessary in order to address 
more fully the time-evolution of the 100 kyr cycle. For instance, the correlation of 
6 18 0 signals from Sites 677 and 607 is intriguing, but comprehensive tests for the 
global coherence of £ 18 0 signals should use data from more than two sites. The 
phase information contained in the paleoclimate records may reaffirm the notion that 
external earth-orbital forcing could be the pacemaker of the ice age cycles. How- 
ever, the mid-Pleistocene amplitude increase of the quasi-periodic 100-ky <5 18 0 signal 
awaits a complete explanation, which appears to depend on factors other than orbital 
insolation changes. 

Future Tasks: 

1) Collection of climate proxy data series, with good time control, from both 
deep-sea, lake-bed, and land-based sedimentary sequences. Much data probably lies 
dormant in the older cores stored bv the Ocean Drilling Program, but only the older 
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cores often suffer from coring gaps that impede detailed time series analysis. Newer 
core data are spliced from parallel- paired drill cores at a single location. 

2) careful cross-spectral analysis of different data series to investigate how the 
earth climate system, as a whole, responded to orbital insolation cycles over the past 
2-3 My, the period of Earth history most relevant to current global change problems. 
The techniques developed for the above project can be extended to cross-series studies. 
The buzzword for this kind of study is ’global teleconnections, ’ the manner in which 
climate at different parts of the globe interacts. The response of climate to modest 
changes in its boundary conditions can (in principle) be tracked by its response over 
individual Milankovitch cycles. 

3) similar analysis of data series from earlier periods of Earth history e.g. the 
Cretaceous and the Eocene. Data is accumulating for a significant response to orbital 
cycles, but a global synthesis has only been attempted so far with numerical climate 
models. The Earth’s climate was warmer in these two periods, but our understanding 
of how the Earth sustained such temperature is, at best, incomplete. Many studies 
point towards higher carbon dioxide levels as a causative factor, and the ubiquitous 
Milankovich-driven limestone/black-shale sequences suggest big changes in ocean cir- 
culation. Are the (cetaceous and Eocene global climates a picture of what awaits us 
in a greenhouse future? 

4) Despite a decade of global- wanning predictions based on numerical climate 
models, it must be admitted that these models (global circulation, energy balance, 
etc) can represent the earth’s climate dynamics only in a crude manner. The features 
of the climate system that are critical to its longer-term variation (decades and cen- 
turies) may not be apparent with the current generation of numerical climate models. 
Improving the global circulation models is a high priority. This includes atmospheric, 
ocean and coupled atmosphere-ocean models. 

5) Collection of a long-term global climate database to ’validate’ the output of 
global circulation models. Such data are essential to check if the numerical climate 
models are working. Comparisons of atmospheric GCM results with satellite data 
have begun, but I am not aware of a global data- validation of oceanic GCMs. 
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Figure Captions 


Figure 1 . Envelope inversions for the time-evolution of the 41-ky signal in the 
benthic <r 8 0 senes from DSDP Site 607 and ODP Site 677, plotted against similar 
analyses tor the astronomical insolation series derived from Berger [19781 (‘65° old’l 
and Berger and Loutre [1991] (‘65°N new’). Seven 47r-prolate Slepian eigentapers 
are used to constrain the estimates, using 'smoothness’ penalty function (10). The 
amplitude units in this and succeeding plots are parts-per-thousand h 18 0 Note the 
phase corrdation j with a constant shift) between the older astronomical series and 
data from ESDI Site 607, and the correlation between the newer astronomical series 
and data from OD1 Site 677. This arises from the orbital tuning performed on the 
data series. Orbital tuning does not enforce correlations between envelope amplitudes 
as is evident lrom the upper panel. ’ 

Figure 2. Similar to Figure 1, but with the envelope of the Berger [1978] se- 
nes omitted and benthic <5 18 0 data from DSDP Site 607 expressed in terms of the 
S ackleton et al [1990] time scale. Seven 47r-prolate Slepian eigentapers are used to 
constrain the estimates, using a ’smoothness’ penalty function. Note the phase corre- 
lation (with a constant shift) between the data series. Note the improved correlations 
between envelope amplitudes in the upper panel. 

Figure 3. lest for the abrupt onset of the 100-ky signal. Double-line envelope 
inversion for the benthic <5 18 0 series from DSDP Site 607 on its published time scale 
at the two major eccentricity quasi-periods of 95 and 124 ky. Seven 47r-prolate Slepian 
eigentapers are used to constrain the envelopes, which are constrained to be smooth 
everywhere except at t a =0.85 Ma. 
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Envelopes for Obliquity Cycle 
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Obliquity Envelopes on ODP 677 Time Scale 
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DSDP 607 - Benthic Oxygen- 18 


Figure 3 
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A First-Order Global Model of Late Cenozoic Climatic Change: 

Orbital Forcing as a "Pacemaker" of the Ice Ages 

Barry Saltzman 

Department of Geology and Geophysics, Yale University 
New Haven, CT 06511 

The development of a theory of the evolution of the climate of the earth over millions 
of years can be subdivided into three fundamental, nested, problems: 

(I) to establish by equilibrium climate models (e.g., general circulation models) the 
diagnostic relations, valid at any time, between the fast-response climate variables (i.e., the 
"weather statistics") and both the prescribed external radiative forcing and the prescribed 
distribution of the slow response variables (e.g. , the ice sheets and shelves, the deep ocean state, 
and the atmospheric C0 2 concentration), 

(II) to construct, by an essentially inductive process, a model of the time-dependent 
evolution of the slow-response climatic variables over time scales longer than the damping times 
of these variables but shorter than the time scale of tectonic changes in the boundary conditions 
(e.g, altered geography and elevation of the continents, slow outgassing and weathering) and 
ultra-slow astronomical changes such as in the solar radiative output, and 

(III) to determine the nature of these ultra-slow processes and their effects on the 
evolution of the equilibrium state of the climatic system about which the above time-dependent 
variations occur. 

In this discussion we touch upon all three problems in the context of the theory of the 
Quaternary climate, which will be incomplete unless it is embedded in a more general theory 
for the fuller Cenozoic that can accommodate the onset of the ice-age fluctuations. We construct 
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a simple mathematical model for the Late Cenozoic climatic changes based on the hypothesis that 
forced and free variations of the concentration of atmospheric greenhouse gases (notably C0 2 ), 
coupled with changes in the deep ocean state and ice mass, under the additional "pacemaking" 
influence of earth-orbital forcing, are primary determinants of the climatic state over this period. 
Our goal is to illustrate how a single model governing both very long term variations and higher 
frequency oscillatory variations in the Pleistocene can be formulated with relatively few 
adjustable parameters. Although the details of these models are speculative, and other factors 
neglected here are undoubtedly of importance, it is hoped that the "dynamical systems" 
formalism described can provide a basis for developing a comprehensive theory and 
systematically extending and improving it. 

The equations for the variations of global ice mass (I), carbon dioxide (ji ) and ocean 
temperature (0), as presented by Saltzman and Maasch (1991) for our model system are: 

^ - «, - a 2 tanh(<» “ V - ct 2 k 9 6 - a^[/?(r) -/?*] + W, (4) 

^ - fix - V + /V - /V - W * + W (5) 

dfi 

-j t - 7, ' It 1 ~ 7 3 0 + W + W 9 (6) 

where c, ke, and k R are constants determined from equilibrium climate model (e.g., GCM) 
experiments relating summer surface temperature at high latitudes to atmospheric C0 2 content, 
to deep ocean temperature, and to the departure of incoming solar radiation at high latitudes, 
R(t), from the present value R*. a 3 and y 3 are inverse time constants for the response of glacial 
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ice and deep ocean temperature assumed to be (lOky)' 1 and (4 ky) respectively, a,, ft, and 
7l are the rates at which global ice mass, C0 2 , and mean ocean temperature would tend to 
increase, respectively, if there were no C0 2 in the atmosphere (n = O), no ice on the planet (I 
= O), no random forcing (W = O), mean ocean temperature was at 0°C, and R at its present 
value R*; these coefficients determine the equilibrium values of I, n, and 0 for any level of 
forcing F and are to be evaluated from the observed late Pleistocene state as an initial condition 
(in a hindcast sense). The remaining six coefficients a 2 , ft, ft, ft. ft. ^d y 2 are considered 
to be the adjustable parameters of our model that can be tuned to account for as much of the 
observed variability and covariability of I, ix, and 0 as possible. 

In Fig 1 we depict in a highly simplified schematic form the interactions implied by our 
model between the three slow-response variables, I, m, and 0, and between these variables and 
both the fast response climatic variables, e.g., surface temperature r, and external forcing due 
to both insolation changes and tectonic variations. A heavy dashed arrow denotes an essentially 
simultaneous, quasi-static (or equilibrated) response of the sign indicated, while a heavy solid 
arrow denotes an inertial time lag in the response. Because cryospheric bedrock depression, D, 
forms the basis of many other models we also include possible interactions of ice load with this 
variable in this diagram. The influence of the slow-response variables on the sign and 
magnitude of fast response climate can be estimated by GCM sensitivity studies and many useful 
results have been obtained; the sign shown on the heavy dashed arrows refers to the influence 
of the slow variables on one particular fast variable, surface temperature r. Although it is 
difficult to calculate the relevant fluxes that determine the slow-response changes, we indicate 
by the signs on the heavy arrows our assumption regarding the signs of the first order effects. 
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In the lower part of the figure we depict by the thinner lines the Bemer-Lasaga-Garrels ( 1983 )- 

type model for the equilibrated response of atmospheric C0 2 to fundamental tectonic and 
weathering processes. 

By assuming plausible time constants for the glacial ice mass and global mean ocean 
temperature, and setting the values of six adjustable parameters (rate constants), a solution for 
the last 5 My is obtained displaying many of the features observed over this time period 
including the transition to the near- 100 ky major ice age oscillations of the late Pleistocene (see 
Fig 2). In obtaining this solution it is also assumed that variations in tectonic forcing lead to a 
reduction of the equilibrium C0 2 concentration (perhaps due to increased weathering of rapidly 
uplifted mountain ranges over this period). As a consequence of this C0 2 reduction the model 
dynamical system can become unstable, bifurcating to a free oscillatory ice-age regime that is 
under the "pacemaker" influence of earth-orbital (Milankovitch) forcing. 

We view this model as an illustration of the potential of a "dynamical systems" approach 

to the formulation of a theory of long term climatic change occurring under the constraints of 
prescribed radiative and tectonic forcing. 
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Radiative Forcing (diurnal, seasonal, Milankovitch) 



Fig. l 


Schematic diagram showing the interactions implied by the model between the three 
slow-response variables /, /z, and 6, the fast-response "climate" variables, such as r 
(inner circle), and external forcing due to insolation changes (upper box) and 
tectonic variations (lower box). The effects of bedrock depression, not included in 

this model, are represented by D. 
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The Orbital Record in Stratigraphy 

Alfred G. Fischer 
Department of Geological Sciences 
University of Southern California 
Los Angeles, CA 90089-0740 


Abstract 

Orbital signals are being discovered in pre-Pleistocene sediments. Due 
to their hierarchical nature these cycle patterns are complex, and the 
imprecision of geochronology generally makes the assignment of strati- 
graphic cycles to specific orbital cycles uncertain, but in sequences such 
as the limnic Newark Group under study by Olsen and pelagic Cretaceous 
sequence worked on by our Italo-American group the relative frequencies 
yield a definitive match to the Milankovitch hierarchy.Due to the multi- 
ple ways in which climate impinges on depositional systems, the orbital 
signals are recorded in a multiplicity of parameters, and affect different 
sedimentary facies in different ways. In platform carbonates, for exam- 
ple, the chief effect is via sea-level variations (possibly tied to fluctuating 
ice volume), resulting in cycles of emergence and submergence. In limnic 
systems it finds its most dramatic expression in alternations of lake and 
playa conditions. Biogenic pelagic oozes such as chalks and the lime- 
stones derived from them display variations in the carbonate supplied by 
planktonic organisms such as coccolithophores and foraminifera, and also 
record variations in the aeration of bottom waters. Whereas early stud- 
ies of stratigraphic cyclicity relied mainly on bedding variations visible 
in the field, present studies are supplementing these with instrumental 
scans of geochemical, paleontological and geophysical parameters which 
yield quantitative curves amenable to time- series analysis; such analysis 
is, however, limited by problems of distorted time-scales. My own work 
has been largely concentrated on pelagic systems. In these, the sensitivity 
of pelagic organisms to climatic-oceanic changes, combined with the sen- 
sitivity of bottom life to changes in oxygen availability (commonly much 
more restricted in the Past than now) has left cyclic patterns related to 
orbital forcing. These systems are further attractive because (!) they 
tend to offer depositional continuity, and (2) presence of abundant mi- 
crofossils yields close ties to geochronology. A tantalizing possibility that 
stratigraphy may yield a record of orbital signals unrelated to climate has 
turned up in magnetic studies of our Cretaceous core. Magnetic secular 
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variations here carry a strong 39 ka periodicity, corresponding to the the- 
oretical obliquity period of that time - Does the obliquity cycle perhaps 
have some direct influence on the magnetic field? 

I consider the following lines of research to be particularly important: 

• (1) St udies of stratigraphic sequences in which Milankovitch cyclicity 
is particularly apparent, and in which the record appears to be un- 
broken and extends for time spans in the 10 6 - 10 7 Ma range. This 
includes such sequences as the Newark Series (Olsen), the pelagic 
Cretaceous of Italy which we have been studying, and the Eocene of 
Angola (Fig. 5). 

• (2) Extending such studies to the tracking of magnetic secular vari- 
ation, which may turn out to provide a record of orbital variations 
independent of climate. 

• (3) Exploring the geographic dimension, by global mapping of the 
distribution of cycle styles for given time-slices. How do cycle pat- 
terns change with latitude, from hemisphere to hemisphere, from 
ocean to ocean? Only such studies will bring cyclostratigraphic 
studies to bear on the problems of climatic change. ALBICORE 
is a start in that direction. 

• ( 4 ) Extending cyclostratigraphic research into the Paleozoic. Mi- 
lankovitch patterns, in particular the 1;5;20 ratio of precession to 
the eccentricity cycles, have now been established back to the Tri- 
assic Period, but Paleozoic stratigraphic patterns do not seem to fit 
this scheme. Were the orbital periods, or the Earth’s response to 
them, different in Paleozoic times? 


1 Orbital Variations 

Quasi-rhythmical orbital variations have affected the Earth since its inception. 
The current patterns of such rhythms with periods of up to 400,000 years, are 
well defined from astronomical observations. Not so clear is how the major or- 
bital parameters and their minor variations have changed through time. The 
length of the day, for example, is lengthening with transfer of angular momen- 
tum to the moon, and the current rate of change has been well established, but 
it seems highly unlikely that the change has been linear, and the existing data 
on this from historical geology are unsatisfactory. Other orbital variations such 
as the obliquity cycle and the precession are linked to the rotation rate, so that 
they too have changed with time, in ways that remain undefined. Astronomers 
are interested in the patterns of change for obvious reasons, but so are geologists 
and climatologists. If the orbital variations of the past have left a record in the 
rocks - specifically, in the sequentially accumulated layers of sedimentary and 
\olcanic rocks that form an incomplete envelope of the crust - they may provide 
a geochronology (Gilbert, 1895) and a means for refining the crude time scale 
provided by radiochemistry. But furthermore, the orbital variations influence 
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the latitudinal and seasonal distribution of insolation, and thereby atmospheric 
climate and oceanic dynamics, and thus come to be agents in ’’Global Change.” 
While the major Icehouse and Greenhouse modes of the outer Earth have prob- 
ably been driven by internal cycles (mantle convection - Wilson cycle of plate 
tectonics: Fischer, 1984). The orbital variations have modulated oceanic an 
climatic behavior within these major modes. Such modulations may be thought 
of as experiments, and if they can be reconstructed from the historical record 
they will bring light to the range of climatic-oceanic behavior that lies beyond 
the realm of human experience. 

2 The Quaternary Record 

The case for such orbital forcing has now been compellingly made for the Pleis- 
tocene. It was first suggested nearly 150 years ago by Adhemar, and the theory 
was further developed and improved by Croll (for a summary, see Imbrie & Im- 
brie, 1979), but its quantitative footing - that the orbital variations vary inso- 
lation substantially - was the life-work of M. Milankovitch (1941), subsequently 
improved by Berger (1980, 1988) and others. The tie of the glacial record to or- 
bital variations did not become definitive until Imbrie and others discovered that 
the isotopic record in the foraminifera of Pleistocene stratigraphic sequences re- 
trieved from the ocean floor provided a proxy of ice volume, and found that 
the fluctuations in ice volume not only showed the same hierarchical frequencies 
of the orbital variations, but also historical coherence between these different 
phenomena (Imbrie, 1982). 

3 Pre-Quaternary Record of Orbital Variations 

But many climatologists and geologists remained dubious about the existence 
of an orbital record during non-glacial times. Duff, Hallam and Walton (1967) 
suggested that whereas the relatively small changes in insolation values during 
glacial times became greatly amplified by the positive feedback of a greatly in- 
creased Earth albedo due to the spread of glaciers and pack ice, the absence 
of such feed-back during non-glacial times made a record of orbital variations 
unlikely. Stratigraphers working in the gap- riddled record of the epicontinen- 
tal regions saw little hope of recovering a record of persistent rhythmicity.Yet, 
some stratigraphers found rhythmic patterns in the stratigraphic record that 
seemed best explained as products of orbital forcing. Thus G.K. Gilbert (1895) 
interpreted the rhythmic spacing of limestones in the Late Cretaceous of Col- 
orado as the expression of the precessional cycle, W. Bradley (1929) viewed the 
rhythmical alternations of oil shales and dolomite beds in the lacustrine Green 
River Formation (Eocene) of the Rocky Mountain region in the same manner, 
and W. Schwarzacher (1947) viewed the alternations of massive and laminated 
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carbonates in the Late Triassic Dachstein platform of the Alps as p recessional, 
and their grouping into bundles of 5 as an expression of the 100 ka eccentricity 
cycle. As geological attention has shifted from purely local or regional con- 
cerns to global patterns, the number of these stratigraphers has grown (e.g., 
ROCC group, 1986; Fischer, 1986; Fischer et al., 1990; Fischer, 1991; Fischer 
and Bottjer, 1991). 

4 Oscillations Recorded in Older Stratigraphy 

The pragmatic facts are that the stratigraphic record is replete with repetitive 
features - some visible to the eye, (Fig. 5), others (such as the Pleistocene 
isotope curve) only retrievable by instrumental studies. Some reflect only the 
stochastically recurring alternations between the several modes of a depositional 
system, such as the alternation of channel and overbank deposits in alluvial 
systems, and were designated as ”autocyclic .” But others seem to have been 
’’allocyclic,” driven by forces outside the regional setting, and candidates for the 
rhythmic climatic- oceanographic changes to be expected from global forcing. 
These oscillations are of many sorts, of which the following have been recognized 
to date: 

• 1. Cryogenic cycles. Changes in global ice volume, reflected in 

- (a) variations in the isotopic composition of sea water. Best recorded 
in foraminiferal tests of pelagic sediments retrieved from the deep-sea 
floor (isotopic cycles) (Imbrie, 1982) 

- (b) oscillations in sea level, on the scale of 10“ C 10 2 m, best recorded 
in subtidal-intertidal alternations and emergence cycles of carbonate 
platforms (emergence cycles), (e.g. Schwarzacher 1947, Fischer 1964, 
Goldhammer et al. 1987, Hinnov and Goldhammer 1991), (Fig. 1). 

• 2. Carbonate production cycles. Oscillations in productivity of pelagic 
carbonate producing organisms (mainly coccolithophorids) are best recorded 
in pelagic chalk and marl sequences (Herbert and Fischer, 1987; Herbert 
and d’Hondt, 1990; Fischer et al. 1991) (Figs. 2, 5). 

• 3. Dilution cycles. Oscillations in the flux of detrital mud are best recorded 
in hemipelagic sediments of continental margins (Roof et al., 1991). 

• 4. Dissolution cycles. Oscillations in the depth and intensity of the lyso- 
cline - the level at which oceanic carbonate accumulation gives way to 
carbonate dissolution, best recorded in relatively deep (2-5 km) pelagic 
sequences. 

• 5. Desiccation cycles. Oscillations in the regional precipitation-evaporation 
ratio are best recorded in 
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- (a) marginal marine evaporite sequences, where annual varving per- 
mits an approach to net evaporation as recorded in annual sulfate 
precipitation (Anderson, 1982, 1984). 

- (b) alternations of lake and playa conditions in lacustrine systems 
(Fig. 1) (Olsen 1987, Fischer and Roberts 1991). 

• 6. Redox cycles. Oscillations in the aeration of bottom waters best 
recorded in pelagic systems by 

- (a) retention of organic carbon (Figs. 2, 5), and 

- (b) shifts in the spectrum of bottom-dwelling animals, best reflected 
in their burrowing patterns (ichnofabric) (Fischer et al., 1991). 

• 7. Magnetic cycles. Oscillations in magnetic parameters may be signif- 
icant in sediments which acquired a remnant magnetism during or soon 
after deposition, and in which this signature has not been irretrievably 
lost by subsequent/magnetic overprints. The remnant magnetism thus 
developed depends (a) on the presence and character of suitable magnetic 
carrier phases (such as the mineral magnetite), and (b) on the strength 
and direction of the then-prevailing magnetic field. 

Inasmuch as the carrier phase is linked to lithology, which responds to climate 
and oceanic change, oscillations in the carrier phase may be expected to reflect 
orbital (as well as other) sorts of lithic forcing. Hence it is not surprising that 
the detailed magnetic investigations of the Piobbico core (Napoleone et al., 
1991, 1992) find the 100 ka eccentricity cycle, dominant in Fourier spectra of 
lithic variation, to be present in the magnetic intensity spectrum as well (Fig. 
4). It is not so easy to explain why it should also appear in the inclination 
and declination spectra. It is even more difficult to understand why a 39 ka 
periodicity - that of the obliquity cycle - should dominate the magnetic intensity 
spectra and should also appear in the inclination and declination spectra, when 
it appears as only a very weak component of the various spectra related to 
lithology. There would thus appear to be a possibility that the magnetic field is 
affected by orbital variations - a suggestion that has been made previously, but 
has never been taken very seriously by the paleomagnetic community. If it were 
to be true, then paleomagnetic studies might provide a record of orbital cycles 
that is independent of transmission through climatic and oceanic dynamics - a 
possibility worth pursuing. 

Paleontological criteria play a large role in the recognition of these cycles (la, 
lb, 2, 4, 6) - which should not be surprising, considering the great sensitivity of 
organic communities to climatic and oceanographic change. 
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5 Theoretical Considerations 

Complications arise from the following factors: 

• 1. The cycles may be overprinted and swamped out by grosser lithic 
changes in response to tectonic-geomorphic events. This is particularly 
the case in marine settings near major sources of detrital sediments, and 
is minimized in carbonate platforms and pelagic settings. 

• 2. Cycles may be only partially preserved or totally lost owing to interrup- 
tions in deposition and continual reshuffling of sediments such as occurs 
in the r tempestite regime.” This implies that many stratigraphic facies 
are never likely to lend themselves to the establishment of a ’’cyclostratie- 
raphy.” 

• 3. Cycles of the higher frequencies may be largely or entirely destroyed 
by the burrowing activities of organisms (bioturbation). This is likely to 
be the case in slowly deposited facies, such as the ’Ted clay” of the very 
deep ocean floor, accumulated at mean rates of lmm/10 3 = lm/10 6 years. 

• 4. Cycle patterns are hierarchical and therefore complex (Figs. 2, 4, 
5). The earliest workers sought to identify stratigraphic cycles with only 
one forcing period, such as that of the precession. Subsequent studies 
such as those by Schwarzacher (1947), Van Houten (1964), Herbert and 
Fischer (198/) found hierarchical patterns. The hierarchy most commonlv 
encountered is the grouping of ca 5 bedding couplets into sets (bundles, 
Figs. 2, 5), which may in turn be grouped into superbundles of 4 (Fig. 2). 
On the other hand, the patterns can become complicated when members 
of the hierarchy shift phase relative to the others, and vary in strength of 
expression. 

• 5. The different orbital forcing functions affect climate and oceanic be- 
havior in quite different ways, and impinge upon a specific depositional 
setting via different pathways. The northern and southern hemispheres, 
for example, respond to the obliquity cycle in phase, but to the preces- 
sional cycle 180° out of phase. When this complication is combined with 
the observation that the climatic-oceanographic forcing of any given de- 
positional setting contains both globally averaged effects such as sea level 
and locally imposed effects such as variations in the amount and tim- 
ing of insolation, the likelihood of a wide range of possible combinations 
and variations emerges. When these effects take different pathways that 
impose different lag times (such as global oceanic turnover), further com- 
plications may result. On the one hand this may be daunting for a first 
recognition of cycle patterns, but on the other such complexities, once 
resolved, provide a wealth of additional information. 
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6 Identification of Specific Cycles 

Vital to development of a cyclostratigraphy is the identification of oscillations 
observed in the stratigraphic record with specific cyclic forcing functions. This 
revolves largely around timing the cycle period. The following approaches have 
been used. 

6.1 Varving. 

The varve method, employed by Bradley (1929), Fischer and Roberts (1991), 
Ripepe, Roberts and Fischer (1991), and Anderson (1982, 1984). Some sedi- 
mentary sequences - in particular those of deep-water evaporites and those of 
meromictic lakes - retain a fine lamination which can with reasonable probabil- 
ity be assigned to the annual cycle, (varving). Continuous varving permitted 
Anderson to plot variations in sulfate precipitation for a 200,000 year record, 
which provided a remarkable record of the precession in late Permian time. 
Episodic varving in lake sediments, extrapolated to the non-varved intervals, 
permitted Bradley to recognize the precession in lacustrine Eocene sediments 
of the Rocky Mountain region (see also Fischer and Roberts, 1991; Ripepe, 
Roberts and Fischer, 1991). Varved sediments are, however, rare, and gener- 
ally do not form time-series long enough to be useful in timing cycles in the 
Milankovitch frequency band. 

6.2 Radiometry 

Radiometric approaches are fairly accurate in the radiocarbon range (the last 30 
ka, possibly extendable to 100 ka), and are applicable to many sediments, but 
for the vast bulk of geological time (Harland et al., 1990) radiometry depends 
on the dating of specific geological events, such as the emplacement of an ash 
layer or an intrusion, which are then extrapolated to the stratigraphy at large. 
Stratigraphical stage-boundaries dated in this manner generally have confidence 
limits of one or two million years for the last 100 Ma or so, but beyond this the 
uncertainties increase toward the 10 Ma level, and, in Cambrian time, beyond 
that. The durations of Mesozoic stages, averaging 3-10 Ma long, have errors 
in the range of 1-5 Ma. Rhythmic time series studied to date generally occupy 
fractions of such a stage, and extrapolating the assumed stage duration down 
to the level of the time-series in question involves further errors depending on 
accuracy of stratigraphic correlations and uniformity of sedimentation rate. As 
a result, such calculations are approximations with confidence limits in the range 
of a factor of 1.5 to 2. This generally serves to ascertain whether a given rhythm 
falls within the confines of the Milankovitch frequency band, but generally does 
not identify it definitively with one of the specific orbital variations. 
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6.3 Magnetic reversal stratigraphy 

Magnetic reversals since the Late Jurassic have now been identified on the sea- 
floor, and can be recognized in many stratigraphic sequences. Through bios- 
tratigraphy these reversals have been tied to the radiometric scale. The width of 
the corresponding magnetic anomalies on the deep- sea floor provides a means 
of refining the radiometric time- scale, assuming relatively constant sea-floor 
spreading rates. At times of frequent reversals, the polarity zones are only a 
few million years long, commonly shorter than stages, and may thus afford a 
better basis for estimating the periods of cycles. Whereas much of our work 
has been in the ’’Cretaceous long normal” polarity chron which lacks the requi- 
site reversals, work in the Tertiary (Schwarzacher, 1987; Herbert and d’Hondt, 
1990) have used reversal stratigraphy to good effect in dating cycles. 

6.4 Ratios 

As pointed out above, stratigraphic cycles commonly occur in hierarchies. A 
grouping of ca 5 bedding couplets into bundles, in Triassic platform emergence 
cycles has now been well documented (Schwarzacher, 1947; Goldhammer et al., 
1987; Hinnov and Goldhammer, 1991). It has been found in Triassic- Jurassic 
lacustrine sequences (Van Houten, 1964; Olsen, 1986), and occurs in Eocene 
(Fig. 5) and Cretaceous (Fig. 1, 2) pelagic sequences (Fischer et al., 1990). 
Furthermore, the Triassic-Jurassic lacustrine beds and the Cretaceous pelagic 
sequence of the Scisti a Fucoidi show a grouping of the 100 ka bundles into 400 
ka superbundles. The geochronology based on radiometric data shows that these 
examples all lie within the Milankovitch frequency band, and thus the case of 
identifying bedding couplets with the ca, 20 ka precession, the bundles with the 
ca 100 ka eccentricity cycle, and the superbundles with the ca 400 ka eccentricity 
cycle becomes compelling. The ratios between cycle levels in the hierarchy thus 
emerge as an important clue to cycle identity. It is noteworthy, however, that 
to date no such good ratios have been found in the Paleozoic. Such studies are 
as yet in their infancy, but the 5;1 ratios, commonly visually striking in the 
Cenozoic and Mesozoic, have not emerged (cf. Boardman and Heckel, 1989; 
Goldhammer et al., 1991). 

7 Present Status of Global Cyclostratigraphy 

At this stage, the case for a pre-Quaternary record of orbital variations has 
been made in principle. The main stratigraphic facies showing heirarchical 
periodicities of the orbital variations are: 

• (1) Deep-water evaporites (Permian, Anderson 1982, 1984). Their varv- 
ing, offers the best age control. They suffer from (a) being too short to 
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apprehend the longer cycles, (b) from difficulties in tying them chrono- 
logically to other facies, and (c) in being scarce. Nevertheless, work on 
these sequences should be pursued. In particular, it now becomes essential 
to restudy the Castile sequence by means of instrumental scans. Studies 
should also be carried on to other sequences of this type, such as the 
varved anhydrites of the Zechstein Formation of Germany. 

• (2) Lacustrine facies: Primarily the Triassic- Jurassic Newark Group se- 
quences studied by Van Houten (1964) and Olsen (1986). Lakes as closed 
systems provide continuity of deposits and a record responding mainly and 
sensitively to local/regional climatic change (wet vs. dry). The disadvan- 
tages of lacustrine studies lie mainly in poor ties to the marine record and 
global geochronology. The most significant work being carried on at this 
time is that of Paul Olsen (Lamont). Other large and persistent lake sys- 
tems of this sort include an unstudied Devonian complex in Nova Scotia, 
which would perhaps provide entry to the presently enigmatic Middle and 
Lower Paleozoic. 

• (3) Biogenic pelagic facies such as those explored in the Piobbico Core 
(Fischer et ah, 1991) in deep-sea cores (Herbert and co-workers). The 
not-so-deep pelagic sediments appear to have recorded (a) changes in the 
aeration state of the bottom waters, and (b) carbonate productivity in the 
surface waters. These parameters presumably reflect changes in circula- 
tion and in the general productivity patterns of the oceans, and a combi- 
nation of local and global effects. Deep pelagic facies are complicated by 
the superposition of dissolution events, and by the effects of bioturbation 
on slowly accumulated muds. 

• (4) Carbonate platform facies such as those studied by Schwarzacher 
(1947), Fischer (1964), Goldhammer et ah, (1987), and Hinnov and Gold- 
hammer (1991). Such facies monitor small-scale sea-level fluctuations - a 
globally integrated signal in contrast to lacustrine cycles. Whereas Mi- 
lankovitch cyclicity has been well substantiated, uncertainties about the 
origin of sea-level oscillations pose a problem (I lean toward small-scale 
glacial effects). Also, like the evaporite and lacustrine records the plat- 
form rocks generally lack the means of close correlation into the global 
stratigraphy, based mainly on pelagic fossils. 

8 My Own Researches 

8.1 Piobbico Core 

I have been working primarily with orbital cyclicity in the pelagic facies - and in 
recent years mainly with the Piobbico core, cut by an Italo-American consortium 
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(Premoli Silva-Fischer-Napoleone) in the mid-Cretaceous Scisti a Fucoidi of the 
central Apennines. We have used this core as a means of exploring various 
techniques of extracting continuous time-series data of various parameters from 
rocks. Figs. 2-4 are a summary of the work to date. We are continuing work 
on this core. 

8.2 Eocene of Angola 

Some of my Italian colleagues, working in Angola, have discovered there what 
appears to be a truly extraordinary Milankovitch sequence in Eocene chalks, in 
which the shale-chalk couplets appear even better defined than in the Albian of 
Italy, as is their bundling into sets of ca 5 (Fig. 5). We hope to make a detailed 
photographic record of these exposures in 1992, and to sample the sequence in 
more detail. 

The regional setting of this sequence - between the extremely nutrient- rich 
up welling belt of southwest Africa and the tropical waters of the Gulf of Guinea 
- may well have provided an ideal site for recording lateral displacements of the 
boundary, of the sort that might be driven by orbital cycles. 

The Eocene, like the Mid-Cretaceous, was a time of Greenhouse Climate, and 
this could well turn out to be the most dramatic expression of orbital/Milankovitch 
cyclicity in greenhouse times. We do not presently have support for this study. 
Eocene time contains numerous magnetic reversals, well tied into the planktonic 
fossil record, and if Angolan sequence retains its original remnant magnetism 
then it should be possible to define the cycle periods with a higher degree or 
precision than has theretofore been achieved. 

8.3 Project Albicore 

It. is one thing to establish the effects and a record of orbital forcing in principle, 
in isolated sequences. Such work may indeed help to define the relative changes 
in cycle periods through time. But they will not shed light on geological prob- 
lems by providing refined chronologies, nor will they illuminate the problems of 
ancient climates. For this it will be necessary to study cyclicity globally and 
for restricted time-slices, which will provide a general view of changes in cycle 
patterns as related to latitude, continent-ocean distribution etc. 

Toward this end I hope to generate a global attack on the pelagic facies 
of one time-slice - the Ticinella praeticinesis subzone of the Albian, about 100 
million years ago, at about the peak of the Cretaceous greenhouse. We chose 
this zone because (1), it shows such striking cyclicity in Italy (Fig. 2, 3, 4), 
(2) it is readily recognized by foraminifera and nannofossils, and (3) it stems 
from a time when high sea- levels left a widely distributed record of pelagic 
sediment. We expect to find such sediments in about 15-20 countries. My plan 
is to carry on studies modelled somewhat (with improvements) on those we have 
developed in Italy (Piobbico core), which will produce comparable data. The 
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work in individual countries would be financed and carried out by groups of 
concerned scientists, advised and aided by an international steering committee. 

Our first step in this direction will be an international workshop, organized 
by Premoli Silva, Fischer and Napoleone, to be held on October 4-9, 1992, in 
Perugia, Italy. The reason for choosing that locale is that it lies within striking 
distance of the outcrops in which our model - the cyclicity in the Scista a Fucoidi 
- can be displayed. 

This workshop will be largely combined with another, the APTICORE con- 
ference of Larson and Erba. This will attempt to organize a parallel project to 
focus on the slightly earlier (early Aptian) events - the eruption of enormous 
quantities of basalt in the mid-Pacific region, and the widespread development 
of oil shales the Selli Bed - wdiich is also very well displayed there. A vital part 
of the participants - that of the world’s Mid- Cretaceous stratigraphers - will be 
equally involved in both workshops. 

The aim of the ALBICORE workshop will include (a) alerting the Mid- 
Cretaceous stratigraphers to the opportunities provided by these approaches, 
which provide a focus very different from the conventional one, (b) educat- 
ing them in the general background, in the need for extended interdisciplinary 
approaches, (c) providing a general forum of exchange on these matters, (d) 
organizing some international ” action groups” w’ho w’ould set out to undertake 
such studies at specific locales, and (e) organizing a supporting organization 
that w’ould provide advice, and support such as providing laboratory facilities 
for specific types of analyses. 

We expect 50-75 people for the combined workshops. Larson has asked 
NSF Ocean Sciences for support, through JOIDES, and this may help to cover 
travel expenses for the 15 or so US participants, but is limited to supporting 
US workers. Our dependence on other nations, many of whose scientists do not 
have travel money, makes it imperative to find funds no thus restricted. I hope 
that the NSF Global Change program will allow for this, but we are likely to 
fall short of support for non-US participants. 
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Five stratigraphic sequences showing hierarchical rhythmicily identified with orbital cycles. 
Upper tier comparison of precessional cycles. Note variations in scale. Lower tier eccentricity cycles. A 
Lacustrine facies, Triassic- Jurassic Newark Supergroup, eastern North America. Aa Precessional cycle 
lake- shore mudstones etc., followed by lacustrine fish-bearing shale (commonly black), succeeded by 
playa mudstones a lake-level cycle. Ab Eccentricity cycles. El,2 cycle: fluctuations in carbonate and 
ana lei me content, attributed to degree of basin Hushing by periodic attainment of overflow E3 cycle is a 
modulation in oxidation level, resulting in alternation of drab and red colors. B Pelagic, Albian Scisti a 
Fucoidi, July. Bs Precessional cycle black, more or less laminated shale (anaerobic) succeeded by 
Chondrites maris (dysaerobic) followed by Planolitcs-bcaring limestone (aerobic). Cycle attributed to 
fluctuations in planktonic carbonate productivity linked to degree of bottom aeration. Bb An El,2 bundle 
of precessional cycles expressed in calcjum carbonate values ( mirror plot ) and in occurrence of black 
shales. Piobbico core, Be Instrumental profiles of 8 m (1600 ka) of Piobbico core, showing darkness curve 
{left) and calcium carbonate curve (right). Black shales in center. High-frequency signal is that of 
precessional couplets; these arc grouped into El .2 bundles, and these into E3 superbundles. (After Herbert 
and Fischer 1986). Bd multi u per spectra of darkness curve and calcium carbonate curve, showing the El, 2 
peak (Park and Herbert 1 987). C Hcmipclagic Coniacian-Campanian Niobrara Formation. Colorado, USA. 
Ca Precessional signal : shale, dark, nonbioturbated, anaerobic?, followed by Chondrites chalk (dysaerobic), 
succeeded by Planoiiies- bearing chalk (aerobic). Cb Precessional couplets defined by detailed calcium 
carbonate profiles and organic carbon content. Berthoud No. 1 State. (Pratt cl al. 1990). Cc Entire Niobrara 
Formation, Berthoud No. I Slate. Left calcium carbonate curve, general; right gamma ray log. Precessional 
cycles not resolved. El ,2 cycles, at lower limit of resolution, are grouped into E3 supcrbundles in sets of 
4, and these into a yet longer cycle which may represent E4 or a still longer (1600 ka) cycle. (After Pratt 
et al. 1990). Cd Fort Hays Member. Adobe Oil Sc Gas Co. Johnson Taylor 1 1 -22, mirror plot of resistivity 
latcrolog. Precessional cycles not resolved. El. 2 cycles grouped into E3 “supcrbundles*. (After Laferrierc 
et al. 1987). 0 Platform facies. Laic Triassic Dachstcin Limestone, Northern Alps. Processional signal: 
massive ncritic limestone containing Lrge clams etc. alternates with perilidai algal laminitcs. Evidence 
that a cusiatic oscillation commonly led to full emergence of platform is furnished by occasional presence 
of relic soils (clayey red to crcc mudstones) and common mudstone-filled 
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Fig. 2— Quantitative expression of hierarchical Milankovitch cyclicity patterns in 
the pelagic Mid-Cretaceous (Scisti a Fucoidi) of Italy (Fig. 1). Left curve: varia- 
tions of gray-scale darkness, by microdensitometry of diapositives; right curve: 
calcium-carbonate values. High-frequency dark-light (low-carbonate/high- 
carbonate) bedding couplets (a,b,c,d,e) represent the precessional cycle (produc- 
tivity and redox cycles combined). A baseline variation in carbonate content (and 
thickness of carbonate beds) groups these couplets into sets of ca. 5, representing 
the ca. 100-ka eccentricity cycles (1,2, 3,4). The enveloping trace shows the 
grouping of bundles into sets of 4, representing the ca. 400-ka eccentricity cycle 
(A,B,C,...). From Fischer et al., 1991. 
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Fig. 5--Pelagic chalk-marl sequence in Eocene of Angola. This is the most 
dramatic visual stratigraphic record of Milankovitch cyclicity known to me, but 
has yet to be studied. By analogy with the Scisti a Fucoidi cycles (Figs. 1, 2 and 
3), I would interpret the high-frequency chalk-marl bedding couplets as alterna- 
tions of (a) calcareous plankton blooms combined with bottom aeration, and (b) 
reduction of carbonate productivity combined with bottom stagnation. This 
appears to be a record of the precession. They are bundled in sets of ca. 5 into 
what would have to be the ca. 100-ka eccentricity cycle. The number of bundles 
in the 20-km strip of coastal cliffs is estimated at ca. 100, which implies a 10- 
million-year record of Milankovitch cyclicity. Study is being planned. 
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Evidence of Orbital Forcing in 510 to 530 Million Year Old 
Shallow Marine Cycles, Utah and Western Canada 

Gerard C. Bond Michelle A. Kominz 

Lamont- Doherty Geological Observatory Department of Geology 
Palisades, NY 10964 University of Texas 

Austin, TX 78713 

John Beavan William Devlin 

Lamont-Doherty Geological Observatory Exxon Production Research Co. 
Palisades, NY 10964 Houston, TX 77252 


We have completed spectral analyses of two sequences of shallow marine sedimentary cycles that 
were deposited between 510 and 530 million years ago. One sequence is from Middle Cambrian rocks 
in southern Utah and the other is from Upper Cambrian rocks in the southern Canadian Rockies. In 
spite of the antiquity of these strata, and even though there are differences in the age, location and 
cycle facies between the two sequences, both records have distinct spectral peaks with surprisingly 
similar periodicities. A null model constructed to test for significance of the spectral peaks and 
circularity in the methodology indicates that all but one of the spectral peaks are significant at the 
90% confidence level. When the ratios between the statistically significant peaks are measured, we 
find a consistent relation to orbital forcing; specifically, the spectral peak ratios in both the Utah and 
Canadian examples imply that a significant amount of the variance in the cyclic records is driven 
by the short eccentricity ('109 ky) and by the precessional (*21 ky) components of the Earth’s 
orbital variations. Neither section contains a significant component of variance at the period of the 
obliquity cycle, however. 

In the Utah example, the spectral separation between the short eccentricity and precessional peaks 
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is longer (by about 30%) than expected for the frequencies of the modern orbital parameters. This 
result is intriguing because a shift in the separation of the peaks by this magnitude is consistent 
with the recent estimate by Berger of the changes in the frequencies of orbital parameters for the 
early Paleozoic. 

We were able to extract these results from poorly dated rocks of great age only after applying 
a new method we have developed, which we call the gamma method, for scaling relative time in 
sedimentary strata. The gamma method works only for well developed cyclic sequences. It makes use 
of the fact that if the cycles are periodic, there will be a predictable relation between the duration 
of the cycles and the durations of the facies that compose the cycles. This relation is tested by 
applying least squares procedures to data from sections measured through the cyclic strata. The 
mast important advantage of our method is that it does not require a high-resolution numerical time 
scale; it can be applied to strata of any age, even Precambrian, provided that they have the right 
type of cyclicity. 

The results of our initial experiments suggest to us that our new method has much promise and 
that the record of orbital forcing in sedimentary rocks, particularly in cyclic strata, may be much 
more robust over geologic time than is generally thought. To test this idea, we are applying the same 
methodology to other examples from older and younger parts of the geologic record. One recent 
result of this work is very encouraging evidence of orbital forcing of shallow marine siliciclastic cycles 
of Cretaceous age in part of the Western Interior Basin in Wyoming. If we continue to be successful 
in this ongoing project, we will be in a unique position to address some long-standing questions 
about the stability of the Earth’s orbital system and about long-term changes in the Earth’s orbital 
parameters. 
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The Cyclic Carbonates of the Latemar Massif: Evidence for 
the Orbital Forcing of a Carbonate Platform during the 

Middle Triassic 

Linda Hinnov 

Department of Earth and Planetary Sciences 
The Johns Hopkins University 
Baltimore, Maryland 21218 


The Latemar Massif is an exhumed carbonate platform in the Dolomites of Northern Italy that 
was deposited during the Ladinian stage of the Middle Triassic. The platform interior is comprised 
of hundreds of vertically stacked, meter scale, subtidal limestone/vadose dolostone cap carbonate 
cycles (’’Latemar couplets”) arranged into upward thinning bundles of 5. The individual couplets are 
interpreted as the products of sea level oscillations occurring with a 20,000 year periodicity, the 5.1 
bundling is indicative of a lower order * 100,000 year modulating component. This combination of 
cycling is interpreted as a platform response to precession forcing with the 5:1 bundling an expression 
of modulation by the eccentricity. 

Spectral analysis confirms the significance of the 5:1 bundling, and identifies other significant 
bundling components related to other components of the eccentricity. The results show that the 
Latemar bundling spectrum is complex, with a splitting of components around the 1:% bundling 
frequency, suggestive of the multiple component eccentricity in the 100,000 year range. On the other 
hand, there is no direct evidence for the major 400,000 year eccentricity component, although there 
is a dominant * 700,000 year cycling in the buildup (a 35:1 bundling period). Other higher frequency 
bundling components may be related to obliquity forcing. Simple models are presented to argue for 
and against an astronomical origin of these cyclic platform carbonates. 
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